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ABSTRACT 


The sensory dynamics of muscle spindles in the toad Xenopus 
laevis have been studied using methods of consistent spectral estimation. 
The receptor discharges at both a constant muscle length, and with small 
random variations superimposed on the steady length, were recorded and 
analysed. Spectral estimates involving the discrete discharges were 
obtained by low-pass filtering the process and applying the Fast Fourier 
Transform to the resulting waveform. Estimates of the power spectrum of 
the receptor discharges at a steady length were used to show that the 
discharges could be adequately described as a renewal process. In addition, 
evidence was obtained which suggests that the power spectrum of the 
receptor discharges at a constant muscle length was acting as the carrier 
for the information related to the superimposed length variations. 
Estimates of the cross spectrum between the applied length variations 
and the spindle discharges were used to obtain estimates of the 'best' 
linear transfer function for the dynamic behavior of the receptor. The 
technique used is extremely well suited to the transfer function analysis 
of simultaneous neuronal spike trains. 

The results show that the variability in the intervals between 
the receptor discharges was such that low-noise recovery of low frequency 
dynamic information was possible by low-pass filtering. Also, the receptor 
was shown to be only slightly sensitive to different levels of static 
stretch. However, in the band of frequencies studied, the sensitivity 
of the spindle to superimposed length variations was strongly dependent 


upon the steady length. 
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CHAPTER 1 


INTRODUCTION 


iO) biology, the quantitative analysis of sensory sytems 
usually involves the interpretation of neuronal spike activity resulting 
from the controlled application of an effective physical stimulus. 

This procedure resembles the so-called ‘black box' problem of electrical 
engineering. Here, the investigator tries, while impressing a known 
voltage or current upon an undescribed system, to relate these, through 
a set of parameters, to the observed fluctuations of these variables at 
a set of responding terminals. The parameters may be nonlinear, as for 
example the dynamic resistance of a junction diode, or frequency dependent, 
as is the transconductance of a triode vacuum tube. The biologist may 
also attempt to describe his systems with a set of parameters. 

Usually though, he must first deal with a number of additional problems. 
In the case of sensory systems he must decide not only what constitutes 
an effective and appropriate input, but also, how to characterize and 
control it. He must also decide how the activity in the sensory nerves 
is related if at all to the applied stimulus. The problem of sensory 
coding has in general not been solved and its interpretation, therefore, 
must be clearly defined in any model which is developed to characterize 
the system. 

Sensory systems usually respond to adequate physical stimuli 

with long trains of nerve impulses or action potentials. Each impulse 
is distinguishable from its neighbors only by its location in time. 
The time course of each is usually short when compared to the interpulse 
interval and each has a characteristic 'all-or-none' nature. The 
intervals between impulses usually exhibit the properties of random 
variables and as a result, many ideas originating from the theory of 
stochastic point processes have become associated with the analysis of 
neuronal impulse activity (1-3). 

Although some of the methods reviewed in (1-3) are extremely 


powerful techniques for characterizing either single or simultaneous 
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spike trains, they can become inadequate when either or both of the 
processes are not stationary in time. In many sensory systems the 

result of the application of a stimulus is to introduce nonstationarities 
or trends into a previously stationary process. In such cases a 
representation of the nonstationarity in the frequency domain becomes 

one of the few practical means of characterization. Normalized with 
respect to input these frequency domain representations or transfer 
functions are now common in the biological literature (4-9). 

Procedures that establish the existence or nonexistence of 
correlation between pairs of simultaneous point processes are clearly 
defined (3). However, the investigator of sensory dynamics is usually 
confronted with the situation where one of the two processes is 
continuous. Here, correlation techniques are not as clearly defined. 

As a result, a number of nonequivalent representations of neuronal 
activity have been used (7, 10-12). Of these, the two most often 
considered are those related to impulse frequency: the number of 

spikes in a unit of time; and those related to interval: the time elapsed 
between a pair of adjacent spikes. Recently McKean et al. (13) have 
obtained evidence indicating that impulse frequency is probably bio- 
logically more relevent. They have also shown that both measurements 
are equivalent at very low frequencies with the interval method intro- 
ducing spurious gain and phase fluctuations into the transfer relation 
as the stimulus frequency approaches and passes the mean discharge 
frequency. A smoothed version of either of these measurements can be 
obtained by applying the stimulus repetitively and averaging the results 
from each application (8-9). Smoothing the impulse frequency measure- 
ment this way results in the production of a function that is propor- 
tional. to the probability density of a nerve spike occurring at a 
specific time after the onset of the stimulus. Smoothing the intervals 
results in an estimate for the average interval, or equivalently, the 
average frequency of spike activity following the stimulus onset. 
Matthews and Stein (8) use both of these methods to obtain their transfer 
function. The average frequency, because of its lower variability, 

was used at low frequencies while the probability density method was 


used at higher frequencies. Groen et al. (14) have also used both 


‘ ; 
fr n oo j 
; x te iin bid had BEI i i 
J in 2 i 
j LF a 
iy 
“he " Cty. j yy 
x { ‘. 
iT ey 
in 
wp 
“ 
A , 
ry f 
ae j 
' 
i 
i | 
' 
< 
5 


he ef i! i 
w 
i f 
t? i 
ta % 
{ 
d i 
\ 5 
n ‘ Din * 
' i J 
r 7, 
i 
r 
ae. a Lm yar 
j 
i 7 ~ 
har i 
h ty é j 
Mie Adal 
p Y } Fe oy : un Ti te, > | 


Ly ee Se yee eo Cee ek he OAAGS: ae py # wuss 


ae eee ee ee ae Mare es ee ay a) 4 | ae 


Pee 7 ay ee aw 9 7 ee ey ie nee “A 


ee: r aR 4 dus | (a Ce 
iis i hi ees 
5 bal: oh wi 


i "7 - alia 


methods in their analysis, however at opposite ends of the frequency 
spectrum to that of Matthews and Stein. 

This thesis is a quantitative study of the dynamics of a 
sensory receptor, specifically, the amphibian muscle spindle. The 
dynamics of a single spindle, isolated with only its sensory nerve 
intact, are described, using broad-bandwidth stimuli, in terms of the 
power spectrum of the point sensory discharge. The results show 
directly what information can be recovered from the receptor output by 
various filtering operations. The potential usefulness of this kind of 
analysis has been discussed by Bayly (15) with reference to the usually 
low pass filter characteristic shown by the synaptic response to a 
nerve impulse. That demodulation of neuronal output resembles a low 
pass filter has been demonstrated in the transfer properties of both 
nerve-muscle preparations (13,16), and nerve-nerve preparations (17). 
The recordings of Fatt and Katz (18), obtained near the end plate region 
of a muscle cell responding to single spikes in the motor nerve, show 
a strong resemblance to the impulse response of a low pass filter. 
Indeed, Stevens (19), in his model-oriented review of Synaptic physiology 
has approximated the postynaptic potential as a single exponential. .A 
description of neuronal ‘integration’ in the central nervous system has 
been given by Katz (20), chapter 10. 

The power spectrum of a point process contains information 
related to both the serial dependence between, and the distribution 
of, interspike intervals and is applicable over all frequencies. These 
properties are not shared by any of the measurements discussed in the 
preceding paragraphs. In addition to the power spectrum, estimates of 
various related spectra, such as the coherency function, have shown 
approximately how 'noisy' stimulus recovery would be from this 
particular sensory receptor under linear filter demodulation. Gain and 
phase function estimates have uncovered the best linear model, ina 
least squares sense, for the spindle dynamics. 

In 1931 Matthews (21) showed that in the frog, stretch applied 
to a muscle, constitutes an effective stimulus to the spindle receptors. 
It was also noted that these organs were not only sensitive to changes 


in length but also to the rate of change of length. As a matter of 
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speculation, Matthews suggested that perhaps structures with 
viscous properties underlying the sensory nerve endings could be 
responsible for this rate sensitivity. Some evidence does exist 
that supports this speculation in the mammalian muscle spindle 
(24) but recent observations, relating to the amphibian spindle, 
are completely contradictory to the idea (22-23). 

Although it was originally intended that a study of the 
role of muscle fiber mechanics in the overall spindle response would 
be included in this thesis, the results of Ottoson and Shepherd 
(22-23) have made this unnecessary. There is agreement that the 
gross mechanical properties of the muscle fibers innervated by the 
sensory nerves do not contribute to the dynamic behavior of the 
receptor. There are however purely elastic differences along the 
fiber which are involved in spindle sensitivity. These observations 
are somewhat contrary to those of Ottoson and Shepherd and will be 
published under a separate cover (25). 

The following two chapters are reviews. Chapter 2 is a 
brief account of amphibian muscle spindle physiology along with a 
description of some of the models proposed for its behavior. Chapter 
3 is a review of power spectrum estimation and describes the impact 
of the Fast Fourier Transform (26) on spectral estimation. Chapter 
3 also consolidates some of the more recent ideas concerning spectral 


estimation that are relevant to the methods of this work. 
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CHAPTER 2 


AMPHIBIAN MUSCLE SPINDLE PHYSIOLOGY 


PP AA Morphology 


Muscle, like all living tissue, consists of the orderly 
arrangement of a large number of individual cells. The stfuctune: of 
each of the cells and the way they are collectively organized is 
specialized to perform a specific function. Skeletal muscle cells are 
elongated forms which usually extend from tendon to tendon in the 
muscle. When excited by Fok ab ey: in their motor nerves, they develope 
tension between the tendons which may result, depending on the load, 
in an overall shortening of the cell. Skeletal muscle cells are organized 
PAtoOAmotoOr units; a motor Unit’) is a group of muscle-cells that receive 
their excitation along branches of the same nerve axon. Motor units 
are also specialized in their function: some are suited for slow 
precise movements, others are better suited for quick-withdrawal or 
reflex-type activity. 

The specialized nature of motor units implies that not all 
muscle cells are identical. pendies on single muscle cells have 
revealed that some cells, when stimulated by activity in their motor 
nerves, contract very quickly but are able to maintain tension for 
only short periods. Other cells contract much more slowly but are 
able to maintain this state over much longer periods. 

Amphibian skeletal muscle has been broadly classified into 
two groups, 'fast' and 'slow'. The division has been shown both 
functionally, by motor units (27-28), and structurally (29-30). Recent 
investigations have shown that functional and structural subgroups are 
present in each of the two major groups Giese However, the essential 
distinction between 'fast' and 'slow' cells remains. The membrane of 
the 'fast' muscle cells is capable of actively spreading contractile 
activity along its entire length via an 'all-or-none' action potential. 
This action potential is initiated in the region of contact between the 
motor nerve and the muscle cell by a nevromuscular transmitter (33). By 


contrast, 'slow' muscle cells are unable to produce action potentials 
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and contraction occurs only in regions of local depolarization around 
motor nerve end plates (34). The larger motor nerves innervate the 
"fast' motor units and form well-localized contacts on the muscle cells: 
the smaller motor nerves end in more diffuse end plates on the slightly 
smaller 'slow' muscle cells (35-36). There is evidence that only one 
region of contact exists between a 'fast' muscle cell and an innervating 
motor nerve, while the 'slow' cell may receive end plates from one or 
scveral motor nerves along its entire leneth (35). The ‘fast’ motor 
units respond to a single impulse in their motor nerve with a charac- 
teristic tension twitch. Repetitive stimulation results in the fusion 
of the twitches into a smooth contraction in which the tension 
generated can be graded by the frequency of the stimulation (cf 32). 

On the other hand the 'slow' motor units develop appreciable tension only 
when stimulated repetitively (37). This tension develops relatively 
slowly and can also be graded by altering the stimulus frequency. All 
mammalian skeletal muscle cells appear to be of the 'fast' type in as 
much as all are capable of producing action potentials (24,38-39). 

In addition to motor nerves, skeletal muscle is also inner- 
vated by sensory nerves. In the frog these axons have terminations 
inside spindle-shaped capsules which, over lengths of about 1 mm, 
enclose a number of small muscle cells in a fluid-filled compartment. 

In the extracapsular regions these so-called intrafusal muscle fibers 

also receive terminations from motor nerves. That the axons innervating 
the capsular regions are sensory has been verified both by histological 
examination (40) and by electrophysiological experiment (41). Contractions 
of the intrafusal fibers have been observed, by light microscopy, to occur 
on stimulation of motor nerves (42). The existence of sensory endings 

on extrafusal muscle cells have been reported in some amphibian muscles 
(42), although these are apparently rare. In the mammal, in addition to 
spindles, the existence of sensory receptors in the tendon region of 

most muscles is well established (43). 

According to Gray (35), in the frog toe muscle extensor 
digitorum longus IV, there are two or three spindle systems. Each 


consists of a bundle of intrafusal muscle fibers with two, three or 
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four discrete encapsulated sensory regions, Only one sensory axon 
enters each capsule; after entering the capsule it divides into a 

number of unmyelinated varicose threads which lie along the intrafusal 
fibers over most of the intracapsular region. A parent axon may 
innervate more than one sensory region. In his more detailed study of 
this region Katz(45) describes a series of nerve bulbs of up to 2-3 im 
in diameter connected by cylindrical tubes about .15jm thick, The 
muscle fiber contains numerous sockets in which these bulbs are seated 
in close-fitting contact. The gap between the respective membranes is 
much narrower than that found at the motor end plates. Katz also describes 
the presence of fine filaments bridging the gap between the membranes. 
In a few cases the bulbs were seen not to be in contact with the muscle 
membrane at all, but floating freely in the intracapsular fluid. In 

the center of the sensory region most of the intrafusal fibers seemed to 
lose their characteristic striated structure and instead appeared — 
reticulate. These ‘reticular zones' were approximately 50-100um long 
and devoid of some 85% of the contractile material present in adjacent 
regions. A few of the smaller intrafusal fibers showed little 
dfirerentiation in the reticular zone* “losing Less than 307 of the 
contractile materials. Katz called the regions immediately adjacent to 
the 'reticular zone' the ‘compact zones'. These areas, each about 250ym 
long, were the site of the majority of nerve bulbs, although significant 
numbers were also present on the ‘reticular zone’. 

Gray's observations (35) provided histological evidence that 
the motor terminations on the extracapsular region of the intrafusal 
fibers were characteristic of motor end plates on the extrafusal fibers. 
He showed that the motor innervation of intrafusal fibers was in fact 
branches of the same nerve axon that innervated the extrafusal fibers. 
Branches of motor nerves that innervated intrafusal fibers formed 
the same type of end plate on the intrafusal fibers that they did 
on the extrafusals. Earler, Katz (46) had also obtained evidence that 
intrafusal and extrafusal fibers were co-innervated. He found that, 
under isometric conditions, stimulation of a single, either ‘large‘ or 
"small', motor nerve, increased the afferent discharge from a single 


spindle while at the same time causing a contraction in the muscle. 
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Sometime later, Eyzaguirre (42,47) demonstrated the difference in the 
time course of the afferent discharge resulting from stimulation of the 
individual motor nerves. On beginning repetitive stimulation of the 
larger fibers there was an immediate acceleration of the discharge of 
the sensory ending. Repetitive stimulation of the 'small' motor fiber 
resulted in the frequency of afferent discharges rising very slowly. 
These observations led him to conclude that intrafusal fibers were also 
of the 'fast' and the 'slow' variety. He was however, able to directly 
detect the presence of action potentials in only the intrafusal fibers 
innervated by the 'large' motor nerves. More recently Smith (48-49) 
has shown that both intrafusal muscle fibers are able to produce and 
propagate action potentials. The fibers supplied by small diameter 
motor nerves however, respond to single stimuli with much slower and 
longer lasting contractions. There is evidence indicating that the non- 
reticulating fibers observed by Katz (45) are members of this group of 
slowly contracting fibers. 

The motor innervation of mammalian muscle spindles appears to 
be separate from extrafusal fibers (50). Sensory innervation is more 
complex in that two sensory nerves originate from the capsule of the 


mammalian spindle (51). 


2.2 Physiological Behavior 

The behavior of single amphibian muscle spindles in response 
to stretching the muscle was first studied in detail by B.H.C. Matthews 
(21). In these experiments the muscle was stretched by loading it with 
weights. During application of the load the frequency of sensory 
discharges increased dramatically after which, under maintained load, it 
gradually declined toward a steady level. When small loads were applied 
rapidly, the frequency of discharge adapted much more quickly than when 
the loads were large. This led Matthews to speculate that perhaps the 
rapid adaptation might be due to differences in viscosity between the 
center and the poles of the intrafusal fiber. On the other hand, he 
speculated, the much slower adaptation to larger loads may be dependent 


on the depletion of some substance from the sensory ending. The peak 


--8- 


iL : 
+ ; 
ae 
be ‘ 
i 
Per | i 
is : : : oF 
- i : a rt 
; : ; ‘ 
! 7 ee ere DAT 
Lo ] ‘ tar pi Sy 
oe ed i ie 1h Aso eal 
ty iw yy = ¥ A LF 
? ; " 
ar OF io 5 es 0 ee f ey a 
‘ 4 oy by 1 iy t (< Pie. * iy 
" i Oat Rees SUE ht Hk Dk Te i a 
Reh ee She : bee Pe Bs pf 
Ar i : ' 
( Wie as ae a 
ao ’ 
1D Sat aa a { 
u A * 
es i 
Ghat i { 
‘ fi 
al 
Wh a 
e 
{ 
‘ 
‘ 
\ 
¥ 4 
¢ 
t 
i i 
i; 
R ” 
4 ye ia 
i, Y : ts 
Bi 
‘ 
- 
. ‘ 
A rai 
, ? 
; oo 
, J 
f a: 4 ‘i tr 
f 4 i % tbe 
* ‘ 
~ i « i t \e iy A Lu) CPR is \ < 
fi 7 “a Ae, 
a 
i 
y cL i 


ite ‘pie 4 
hy ’ ; ae i, 


hay 


i bi ee 
‘a 


, aot by o 0 cA 
ag Ho 19.0 i ra 


hey 


/ 7 nn ayy ie 
aH eile ent. Seals he 64 


; 


ce OMT AY 
Salat 
Ne 


frequency reached during application of the load was proportional to the 
rate at which it was applied. The adapted discharge rate was approximately 
proportional to the logarithm of the maintained load. Matthews later 
showed (52) that the afferent discharge ceased when the whole muscle 
contracted, suggesting that the sensory ending was in parallel with the 
fibers that produce the bulk of the tension. Sometime later, Katz (Ab), 
using controlled stretches of constant velocity, showed that the response 
of the spindle to applied stretch depends on two main factors: the 
velocity and the amplitude of the imposed increase in muscle length. 


‘static’ 


He divided the response to these variables into a 'dynamic' and 
component respectively. Recently, Shepherd and Ottoson (53) have more 
quantitatively related the response pattern to the parameters of the 

ramp stretches. 

The variability of amphibian muscle spindle discharges has 
been examined by Buller et al. (54). Their results indicate that at 
very low discharge frequencies the standard deviation in the frequency 
approaches the mean. Hence, in the statistical sense, the discharges 
of the spindle at very short lengths resemble a completely random 
process. As the mean frequency is increased, by stretching the muscle, 
the standard deviation increases to about 3.5 impulses per second ati a mean 
rate of 8 impulses per second. Further increases in mean discharge 
rate do not result in significant changes in the standard deviation. 

The speculation is that this behavior may be due to a constant noise 
voltage due to thermal agitations in the fine sensory terminals. 

Stimulation of the motor nerves which innervate a spindle can 
result, as has already been pointed out, in the acceleration of the 
sensory discharge (42,46,47). Matthews and Westbury (55) have shown that 
repetitive stimulation of 'fast' motor fibers markedly excited the 
spindle when the muscle was at a constant length. The dynamic response 
to stretching however, was left almost unaffected. On the other hand, 
pepe e stimulation of the 'slow' motor fibers only weakly excited 
the sensory ending when the muscle was at a constant length, while 
stretching the muscle resulted in an augmented dynamic response. As 


Matthews and Westbury point out, this is in direct contrast to the effect 
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of the large and small fusimotor fibers on the mammalian spindle. 


2.3 Mechanisms Underlying Receptor Behavior 

In 1950, Katz’ (41) showed that when the recording electrodes 
were placed very close to the spindle capsule, the afferent discharges 
were superimposed upon a local potential which re~developed after each 
discharge. By applying a local anesthetic to the preparation he he 
able to abolish the spike activity and record undistorted the potential 
charges near the sensory ending. The frog spindle thus became the first 
mechanoreceptor from which a receptor potential was recorded. Katz 
observed that the size of the receptor potential varied continuously 
with the amplitude and the rate of stretching. He therefore considered 
the receptor potential to be the immediate cause of the sensory spikes 
present in the unanesthetized condition. He also suggested that the 
conversion from the mechanical stimulus to electrical energy in this 
receptor does not involve an ‘active response' of the nerve endings, at 
least not the kind that is associated with the action potential. 
Sometime later Edwards (56) showed that the potential change produced by 
an external current applied to the sensory axon near the capsule would 
sum with the receptor potential. He was able to show that the discharge 
frequency resulting from stretching the spindle could either be raised 
or lowered depending on the direction of the applied current. This 
was direct evidence supporting the view that the generator potential 
was the underlying cause of the sensory spikes. 

Katz (41) was also able to make measurements relating the 
magnitude of the receptor potential to the frequency of afferent spikes. 
His data showed a linear correlation coefficient of .97 and was thus 
significant evidence that the two were linearly related. Edwards (56) 
using external currents also obtained a linear relationship between the 
applied current and the frequency of discharges. He also noted that the 
spike producing mechanism showed no adaption to step applications of 
current. When the direction of the applied current was reversed, the 
relationship between the current and the discharges showed a markedly 


different slope. This was attributed to the rectification properties 
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of the nerve membrane. Ottoson and Shepherd (57) have shown an 
approximately linear relationship between the following: rate of stretch 
and rate of rise of the receptor potential, rate of stretch and the 
dynamic peak attained by the receptor potential, the amplitude of the 
stretch and the magnitude of the partially adapted receptor potential. 
These results were obtained by applying ramp stretches directly to the 
encapsulated region of the spindle. A linear and instantaneous relation 
between applied current and discharge frequency has also been shown in 
the mammalian muscle spindle (58), although adaptation to constant 
currents seems to be present in some crustacean stretch receptors (6). 
In a paper which preceded his study of the receptor potential, 
Katz (59) described the presence of small 'all-or-nothing' spikes which 
were quite distinct from both the receptor potential and the full sized 
action potentials. These were present in greatest numbers when the 
muscle was under low tension and tended to disappear as the muscle was 
stretched. They were observed to occur in a number of different sizes 
and also formed the prepotentials on the main spikes. Katz suggested 
that these small spikes were formed in the terminal branches of the 
sensory axon and sometimes failed to initiate active propagation past 
the branch points. Those that failed to initiate propagation were 
observed as abortive spikes. The maintained receptor potential at 
eréater muscle tensions probably facilitated the passage of the active 
potentials through the branch points. Recently Ito (60) has produced 
a quantitative report which describes the effect of muscle length upon 
the pattern and frequency of the abortive spikes. The ideas here generally 


support the speculations of Katz(59). 


2.4 Models 


Until recently, most evidence related to the mechanisms of 
muscle spindle dynamics, both in the amphibian and the mammal, pointed 
to muscle fiber mechanics as the main underlying factor. Matthews' (21) 
original suggestion that the form of the afferent discharge in response 
to stretch is a reflection of the differences in viscosity between the 


center and the poles of the intrafusal fiber, has received much 
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speculative attention. In a theoretical treatment Katz (41) related 

the possibility that changes of electrical capacity in the terminal 
nerve membrane, as a result of stretch, could contribute to the phasic 
response of the spindle. He showed however, that this mechanism could 
not provide a phasic receptor potential of the observed magnitude unless 
there was a substantial mechanical amplification of the stretch. In 
introducing his later histological work (45) he said “it is natural 

to suppose that a mechano-receptor such as the spindle owes its high and 
specific sensitivity, at least in part, to the way in which the mechanical 
stimulus is bought to bear on the membrane of the sensory nerve terminal." 
The results of this work showed that most intrafusal fibers showed a 
"reticular zone' in which there were few contractile filaments. 
Excitation of the motor nerves would therefore stretch this zone and thus 
explain the observed increase in afferent discharge (42,46,47). Also, 
the polar region, with its mass of interdigitating contractile filaments, 
will probably be much more viscous than the ‘reticular zone' where these 
filaments are largely absent. Therefore, it might be expected that 
changes in muscle length would first appear across the ‘reticular zones', 
followed by a gradual and partial shift to the polar regions after the 
muscle length is fixed. In addition, the works of Edwards (56) and 
Lippold et al,(58) suggest, because of the failure of the encoding 
mechanism to adapt to constant currents, that mechanical factors are 
very prominent in overall spindle dynamics. Indeed, there is direct 
evidence that mechanical factors are involved in the adaptation of 
mammalian muscle spindles (24), in the mammalian Pacinian corpuscle (61), 
and receptors in the skin of the frog (62). 

In light’ of the evidence, a number of models have been developed 
for both the amphibian (63-66) and the mammalian (66-69) muscle spindle. 
Of these, the model by Houk et al. (65) seems most typical and will be 
briefly described. This model assumes that the 'reticular zone‘ of 
the intrafusal muscle fiber is less viscous but stiffer than the 
adjacent polar regions. Also, the viscosity of the ‘reticular zone’ and 
the stiffness of the polar regions are assumed negligible, The intra- 
fusal muscle fibers are therefore represented as a spring in series with 


the parallel combination of a force generator and viscous dash pot. 
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The spring represents the lumped stiffness of the 'reticular zone' while 
the dash pot and force generator depict the lumped viscous and contractile 
properties of both polar regions. The effect of stretching the model 

is to produce strain in the elements representing each of the regions. 

The resulting generator potential is considered proportional to the 

strain in each of the passive elements weighted by the number of nerve 
terminations present on each of the regions in the actual muscle. Katz 
(45) reported that on average, 3.5 times as many nerve bulbs were located 


" than there were on the ‘reticular zone'. Strain 


on the ‘compact zones 
of the polar regions will therefore be more prominent in the generator 
potential when the model is passively stretched; while strain of the 
"reticular zone' will be the sole contributor to the generator potential 
during motor nerve input under isometric conditions. The encoder in 

the model produces a train of impulses whose frequency is directly 
proportional to the instantaneous magnitude of the generator potential. 
The parameters of this model have been adjusted so that the theoretical 
responses resemble those of the amphibian spindle (66). 

In 1965 Buller (70) published a model designed to test the 
suppositions of the earlier work of Buller et al. (54) concerning the 
events occurring at the sensory nerve terminals of the frog muscle 
spindle. These experiments had shown that, above a certain minimum 
mean frequency, the standard deviation in the frequency about the mean 
was approximately 3.5 impulses per second. Below this critical 
frequency, the standard deviation tended toward the mean as the mean 
decreased. It was suggested that the critical frequency corresponds to 
threshold depolarization of the nerve terminals with further depolariz- 
ation resulting in directly proportional increases in the discharge 
frequency. Variations in the frequency are introduced by a source of 
random potential variations due possibly to molecular agitation in the 
mechanical receptor substance or to ionic noise in the terminal nerve 
membrane. Discharge below the critical frequency would therefore be 
the result of noise transients on the subthreshold depolarization. 
Buller's model was based on these assumptions and its responses show a 


remarkable similarity to the results obtained earlier from the actual 
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spindle. In addition, the results show that the variability which exists 
in the discharges of the model above the critical frequency are not 
particularily sensitive to changes in the rms value of the disturbing 
noise used to represent the nerve terminal fluctuations. Changes in the 
shape and the bandwidth of the noise spectrum have even less effect upon 
the variability. It was therefore concluded, that, in spite of the 
probable changes in the terminal membrane electrical time constant brought 
on by the different levels of depolarization, the disturbing noise may 
be considered as remaining constant over the physiological operating range 
of the receptor. 

More recently, Stein (71) has pointed out that the distribution 
of intervals collected by Buller et al. (54) from the single frog 
spindle is extremely well described by the well known Gamma density 
function. Hagiwara (72) has shown that adjacent intervals from muscle 
spindles in the sartorius muscle of the toad are statistically independent. 
Apparently then, the discharges from the amphibian muscle spindle can 
be considered as realizations of a renewal process and therefore, under 
steady state conditions the Gamma density would be a complete description 
of the process. The Gamma density has also been used to describe 
empirically the maintained discharge originating from the retina of the 
cat (73) and the responses of primary auditory neurons to acoustic 


stimuli (74). 
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CHAPTER 3 
SPECTRAL ESTIMATION 
eel patroduetion 
The power spectrum $4 bof) of a random process as a function 
of time t and frequency f is defined by 


—~j2nf£T 
? dt Gr) 


Cc 
= R Ap Jeet 
where RL 6Es 7) is the autocovariance function of the process. The 
autocovariance function is the mathematical expectation of the product 
of two values which the random process x(t) assumes for two instants of 


rames, Lnat. iss 
Ri 6b 27) = Est oCe)xCtkt)} C22) 


defines the autocovariance between the value of the process at time t 
and its value at another time ttt. For a stationary random process 
RL (tt) can be shown to be (75) 


Aimy ced, 
Too OT Jot x(t)x(t+t) dt (393:) 


R Get) = a tT) = 
eG) a 

The cross spectrum between two processes x(t) and y(t) is defined 
as 


=IoTLe 


Sy (tf) = Jr Ryy (ts Te dt (3.4) 


For a stationary random process Ra Meene is given by 


Pia ee b 

= = Ta 5 SS) 
EB ae He ne ME, Ge ACR a meal a (a) 
S tof) and § Fiona for a stetionary random process, become simply 

< x 

S (f) and S$ (f£) respectively. By manipulation of the expressions 

xe xy 


which define these spectra it is possible to show that, among other 
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things, Se is a real function of frequency while eee is usually 
a complex function of frequency (75-76). 
If the process y(t) is derived from x(t) by passing x(t) 


through a linear system, it can be shown that 


Sees = H(f) S(t) (376) 


and 


(f) 


aA jH(E) |? s(£) (3.7) 


where 


=< 2a t 


H(E) = i h(t) e dt (3.8) 


h(t) being the response of the system to a unit impulse. H(f) is the 


frequency response of the system, 


The squared coherency spectrum between two stationary processes 


(0) and: y Gt yeas: defined as’ 


Is, CO? 
Vee Ee a 79) 
xy S x FD SIGS) 
it'ean be’ shown that 
2 Sy ho) 
US a hice (Ee <a ( ) 


The coherency spectrum is a measure of the linear correlation 
between the processes x(t) and y(t) as a function of frequency. For 
the linear system described by equations 3.6 - 3.8 ee) = AV EOr aliee, 
Bendat and Piersol (76) describe the situation of a squared coherency 
greater than zero but less than one as indicating one or more of three 
possibilities: 

a) extraneous noise is present in the measurements; 


b) the system relating x(t) and y(t) is not linear; 
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c) y(t) is an output due to an input x(t) as well as to other 
inputs. A coherency of zero at all frequencies implies that the 
processes x(t) and y(t) are completely uncorrelated. Goodman (77) 
describes the squared coherency spectrum for a nonlinear system as a 
measure of the degree to which the output y(t) at a particular frequency 
is related by a linear time invariant operator to the input x(t). The 
quantity Spy F)/8, ©) is the linear operator that best approximates 
y(t), in a least squares sense, by acting on x(t). 

The squared coherency spectrum has the interesting property 
that makes it invariant under linear filtering operations. For example, 
Lf Vey Cf) is the squared coherency spectrum between the processes x(t) 
and y(t) and ee) is the same function between u(t) and v(t) then 
a) = yee ) if u(t) is derived from x(t) by passing x(t) through 
some linear filter and v(t) from y(t) by passing y(t) through some other 
linear filter. 

The most serious problems associated with obtaining reliable 
spectral estimates are those which result from the application of the 
previous formulas to finite pieces of data. Today, two practical 
approaches, each with distinct advantages, are available and can result 
in good estimates of the various spectra from finite lengths of data. 
Neither of the methods, however, is entirely clear cut and usually must 
be applied several times with systematic modifications to produce the 
final estimate. A rough idea of the shape of the spectrum before 
estimation is extremely helpful in reducing the labor involved. 

The oldest of the two approaches to power spectrum estimation 
is that formalized by Blackman and Tukey (78). Until 1965 this approach 
was used almost exclusively and offered the additional advantage that 
an estimate of the covariance function was a byproduct of the operation. 
In 1968 Jenkins and Watts (79) extended these ideas to the estimation. 
of the cross spectrum between several processes. The subsequent role 
of the cross spectrum in the estimation of the squared coherency spectrum 
and the gain and the phase functions was also developed. 

In 1965 a highly efficient algorithm for the digital evaluation 


of complex Fourier coefficients was discovered by Cooley and Tukey (26). 
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The Fast Fourier Transform has, in the minds of some (80), revolutionized 
the processes of spectral estimation. The so-called direct route to 
spectral estimates is computationally several orders of magnitude faster 
than the older Blackman-Tukey approach. On the other hand, the direct 
method suffers from the disadvantage that the covariance function is 

not a byproduct of the operation. The role of the Fast Fourier Transform 
will be discussed further in the last section of this chapter after the 
ideas underlying each of the approaches to spectral estimation have been 


described and compared. 


3.2 The Sample Spectrum 


Estimates of the power spectrum of a stochastic process are 


derived from the so-called sample spectrum defined as 
A dae =e 3.11 
Sa a ek eee dt (19) 


where R_(t) is the estimator for the autocovariance function Rt) 
xX 


Two widely used estimators for the autocovariance function are 


ae [a ate 


aOR GOIN iste SE acces cu) dt <|t|<tT 
0) Loin 
Rea) eS (35125 
[Bae i ees pelle cpm, (etl) dee 0</ 7 (<0 
0 [ct |>T 
where 
x(t) O<t<T 
Beis) = (3.13) 
: 0 T<O5 27 


The R(t) given above are each functions of the random variable 
XX 
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x(t). Therefore GD) is also a random variable. This means that 
Be) at all t will be subject to some probability density function. 
This density function is called the sampling distribution of the random 
Varraple Re. C(t). 

XX . ‘ 

It is useful to compare the estimators Roy © and Rx 
by calculating what are usually called the sampling moments of their 
respective distributions. It has been shown (79) that the expected 


values (means) for the estimators are given by 


- Gy fae foler 
B{ R(t) } = (Sen) 
ri lc |>T 
and 
‘ Ree) [z|<t 
E PR) yg (3.15) 
0 fol 


Thus R 2 6) is what is called an unbiased estimator of REO) whereas 
2 x 


R wD is only asymptotically unbiased as the record length T tends to 
Be 
infinity. The variances of each of the estimators are 


Tone yy ees RoR Oe REO fees ia dr 


|r| <2 (3.16) 
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It is seen that these expressions are identical with the exception of 
the multiplying factors preceding the integrals. The variance of 

R 96) grows without bound as |< |>r. The mean square errors in each of 
xX 


ie 


Mis dt nH Ata 


Ani 


the estimators are 


n mae . 2 
E Cie CAR cael Var ee i Rd (34.108) 
E (IR (t)-R,, (0) 1215 van {R. (t)} (3.19) 


A 


It is usually true that the mean square error in Rep is less than 
that, in A al: is ; & ie xe 1 = 
Rog herefore, Roy ho) is preferred to Rg &) for auto 
covariance estimates. In the remainder of this chapter ea = R. 1 
is assumed. 
By using the definition of the sample spectrum given by 


equation 3.11 it is possible to show that 
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= 7 xp) Xe (£) 


le 


x, CE) |? (3.20) 
where XC) is the Fourier transform of Xp Ct) and Xp* CE) is its complex 
conjugate. x (£) can be expressed as the sum of real and imaginary 
components 
= EY j Sw A: 
Xf) Ap) oe jB, Cf) ( ) 


where j=v-1. It can be shown that 


ie x(t) cos 2nft dt @e22) 


ACE) 


and 


B.(£) = [g-xp(t) sin 2nfe de (3.23) 
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Therefore 


ues = + [A,2 CE) pba (224) 
It has been shown (79) that if x, Ct) is a sample from a Normal 

process then, since the Fourier transform is a linear operation, A, (£) 

and B,.(£) are also Normally distributed. Therefore Su Gs is distributed 

as a chi-squared random variable with two degrees of freedom since it is 

the sum of the squares of two Normal random variables. If the sample 

x7 (t) is not Normal, A, (£) and B, Cf) tend to Normality as T becomes 

large as a consequence of the Central Limit Theorem. Hence the 

distribution of $_(E) will be very nearly as chi-squared with two degrees 

of freedom irrespective of the distribution of x(t). 


The expected value of the sample spectrum estimator is 
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Hence, using the convolution theorem 


a 
E {s_()} = es T Ee S_ (f-8) dg 
ae) Sno (oaic (fee jade (3% 26) 


The function W(f) is a slit whose width is of the order of 1/T so that 
for large T it is reasonable to assume that 8 x is approximately 


constant over the slit. Hence 
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sufficiently large. 
As the length of the available record becomes larce, tiican 


be shown (79) that the variance of the spectral estimator tends to 
Var oe MCS) Stal (3.28) 


For a white Normal process this result is exact. Thus the variance of 
Bye) is dominated by a term which remains finite as the length of the 
pecora tends, to infinity. \WjPhis\ indicates that Sn (FD LS NOt 2 Consistent 
estimator of S() in the sense that its distribution does not tend to 
cluster more closely about the true spectrum as the sample size increases. 
It is almost equally probable for ue to lie anywhere from 0 to 

28. fF). The reason for the inconsistency of Sn FD can be described 
simply as follows. The estimate of the spectrum at frequency f is 
actually drawn from a band of frequencies W(f) whose width is about 

1/T. As T is increased the power at f is estimated over narrower and 


narrower frequency bands. However, the efficiency of the estimate of 


the power in the narrowing band does not improve. 


3.3 Consistent Power Spectrum Estimation 


The inconsistency of the sample spectrum has forced investi- 
gators to seek other means of estimating power spectra. The technique 
universally adopted is to trade the 1/T resolution afforded by CE) 
for the lower variances provided by an estimator with a coarser 
resolution. However, not all techniques that result in a loss of the 1/T 
resolution produce a corresponding decrease in variance (aL); 

The Blackman-Tukey or indirect approach to obtaining consistent 
power spectrum estimates is to consider only that part of the auto- 
covariance function to lag t=M where M<<T. This procedure can be thought 
of as looking at the estimator of the autocovariance RS tt through a 
rectangular window of width 2M, or more simply, as a multiplication of 


the autocovariance estimator by a function Wy (tT) where 


-22- 


if 4 | oe AS 
‘ a a 
; A ie 
i os i a 
4) : i 
f Ny ps 
; D 1 eae iy 
ow Ae ; 
, { i vt 7 
1 p: ; ; be 
ps Al ™ “7 
Pal beat ; ) aa 
te i 1 i : 
1 Srey | Plo , 
; p . ; ‘a 
No) ee ala 
i Bi r j ued Mad pif! 
ie i ‘ ’ 
2 s 
iN B Cry; } iy ‘ ; 
i} : } , 
\ 4 1 ‘ mc} 
P| { t . oy % 
i 
* 
nm 
, ' ‘ 
aa ‘ f 
( } Ph ke ‘ ee 
a ae : } 
) 
: 1 1 
re a 
Y 
' F 
2 
’ 
) 
4 
1 
: ey i i ¢ 
y 
' 
4 
F i : 
b 
, 
x 
a 
. “ 
LS ie 
Yr p' ~ 
' &: 2 f 
U ew 4 
“ar ‘ tithe 
£ . 
rs ‘ 
on a a! 
; $0’ b4-d o ? hi ‘ 
v are 
(* af 
a i. 4 yf es ik | 
Ly or iS \, y 
rey ; ae | 6 
* to fat 
a Wh "ep 7 A 
{ Tie ‘poe 
Ma ye . 
t i» - 
"a ks nih F Lon 
ay vat a . a 
I 14 oa ro: ; i i} 2 3 ate, 
? x (tee 
fr ; 
whe 3} } 
rh 
i 
; a eee 
: f a Perri have stn 
4 : ral } ¢ 
yt ~ ie |. ei A) ol enon 
ie © he i SM pa d oh age aif iA Ye ne 
bay : ia y 
bet ars 
bao ‘ bidet (Z : on fds SH. oi i eae i crs ut Vite 


J Sig ‘Ss ithe ; Hae on & i yt 


Reery ane he cs 
Peon Wr a 


ee 


oy j & iD wl 


i iy ial i 7 oh, iy 


0 
ay) ae 
g 


1 |t|<M<<T 
toe Gay (3.29) 
0 | x |>™ 


Blackman and Tukey (78) refer to Wp (t) as a rectangular lag window. 
The smoothed power spectrum estimator is defined, using the 


convolution theorem, to be 


S(f) = [7 Mg (e)8_ (E-g) de (3.30) 


Wp (£) is the Fourier transform of Wp (t) and 


[ sin2n£M 


Watt) = 2M | “Tag 


| —2< F<o G6) 
and is referred to as a spectral window. The expected value of the 
smoothed estimator is (79) 


vat 


Econ) sae ey (=) cutee) ade CLy 


This expression is very similar to that derived for the expected value 
of the sample spectrum with the important difference that Wp (£) now has 
a major base width of 1/M. This means that estimates at a particular 
frequency will be influenced by power at least 1/2M away. On the other 
hand, if SE) is slowly varying in W(t), s__(E) tends to be an 
unbiased estimator of the true spectrum. The variance of S ee (FD ean 

be shown to be (79) 


“aw 


me 2 a 
Var {s_@) Pes cae ie w(t) dt 


a 2 (£) (3.33) 


The variance of the smoothed estimator is reduced by a factor 2M/T 
from that of the sample spectrum. Also, the variance tends to zero as 


T becomes large. 5S mee is therefore a consistent estimator of Sf): 
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The direct approach to consistent power spectrum estimation 
is to break up the data record of length T into subrecords of length 
L with L<<T. The sample spectrum is computed for each segment via the 


relation 


a 2 
s, @) =elx, wl (3.34) 
k k 
where k refers to the kth segment. The smoothed power spectrum 
estimator becomes 
t 
k=L 
Cee ee (3335) 
Sel) T k=l °xx 7 
It can be shown (82) that the expected value of this estimator is 
Ess (iy) al ec We Gayo (ioe Jade (3.36) 
with 
2 
sintLft 
Wf) = L fara (3237) 


aR 


Wp (£) is known as the Bartlett spectral window. The variance of Sf) 
can be derived heuristically. From the previous discussions, it is 
known that the variance of each of the i is approximately equal 


to S2 (£) regardless of the length of the subrecord. Since the eae’, 
XX 


are T/L independent estimates of Sf) it follows that the variance of 
Smuce ls reduced oto, b/T that, Or Sekt). sbateds 
XX xX 
Var fs (Epes dg? (£) - (3:38) 
xX i axx 


It is seen that if L=2M the direct and the indirect approaches both 


result in smoothed consistent estimates with approximately the same 


varlance,. 


The bias in the smoothed estimates can be compared by looking 
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at W(t), the spectral window corresponding to the rectangular lag 
window, and Wp), the Bartlett spectral window. To make the comparison, 
the variances of the two smoothed estimators were set equal, that is 
with L=2M. When using the Blackman-Tukey approach, any bias in sf) 
is due to leakage from neighboring frequencies through Welt); while with 
the directly approach, leakage is through W, (Cf). Half of each of these 
windows, centered at zero frequency, is shown in figure 3.1. Both 
estimators then, result in good estimates of the power in the band 
about some frequency f. Both approaches allow, by a simple choice of 
paOre, Chie’ possibility Ofvexchancsine spectral. resolution for a 
corresponding decrease in variance. As far as bias considerations are 
concerned, if the true spectrum is reasonably smooth there is little 
to choose between the two approaches. If however, SF) contains 
large fluctuations, the large and slowly declining side lobes of Wp (£) 
can result in widespread leakage between regions of different power 
concentration. Large bias in spectral estimates is often not acceptable 
and a process called quadratic modification is used to reduce its 
effects. | 

Quadratic modification is the process of viewing the auto- 
covariance function through lag windows other than the rectangular 
window. The windows are chosen such that their spectral equivalents 
have lower side lobes and die away faster than those of the rectangular 
lag window. There are many lag windows that have these properties 
(78,79). One of the most popular means of reducing side lobe leakage 


is via the Hanning lag window 


1/20 cos ao | «| <M 
w(t) = (3-39) 


O |< |>M 
The Hanning spectral window is shown compared to the spectral equivalent 


of the rectangular lag window in figure 3.2. In addition to increasing 


the rate of side lobe decay, the Hanning window has a considerably wider 


95 


_ FREQUENCY 


Pisure! 3eL 


SPECTRAL EQUIVALENT OF THE RECTANGULAR LAG WINDOW, Wp Cf); 
AND THE BARTLETT SPECTRAL WINDOW, W,(f). 
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Figure 3.2 


SPECTRAL EQUIVALENT OF THE RECTANGULAR LAG WINDOW, W(t), 
AND THE HANNING SPECTRAL WINDOW, Wy). 
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major lobe width. This increased loss in spectral resolution is of 
course accompanied by a decrease in the variance of the estimator. It 


can be shown that (7/9) 


f Sx) 1 
Var 1S GE) e re eG, dt 
(3.40) 
Ay 75M. oo 
i _ Six) 


as opposed to a multiplying factor of 2M/T with ene rectangular lag 
window.’ Reductions in variance from L/T sé Cf) can also be obtained for 
a fixed L and T with the direct approach. This is achieved by averaging 
in the frequency domain (84) or by overlapping the subrecords (83). 

The process of linear modification which has arisen with the 
increasing popularity of the direct approach has been shown to result in 
a loss of spectral resolution without a compensating decrease in variance 
(81). Linear modification is the process of looking at the raw data 
through windows other than the rectangular data window. All lag 
windows can, with minor modifications, be used as data windows. Some 
authors (80,85) feel that some minor linear modification is desirable as 
this process is extremely effective in reducing the leakage of power 
between adjacent frequency bands. 

The distribution of ae) was described earlier as being 
approximately that of chi-squared with two degrees of freedom. The 
two degrees of freedom arise from the condition that a) is the sum 
of the squares of two Normally distributed random variables. Strictly 
speaking, chi-squared with k degrees of freedom is the sum of k Normally 
distributed random variables each with zero means and variances of one. 
Therefore, some normalization of sie) is required since the random 
variables A, (£) and B,C) do not necessarily satisfy these requirements. 
Tt can besshown, ¢/9) ‘that sh ‘the mean of x(t) is removed and if f is 
somewhat removed from zero then, 28 F)/S_ ) is approximately 
distributed as chi-squared with two degrees of freedom. It can also 


be shown that §S ahh) is also a chi-squared random variabie but with 
x 


greater than two degrees of freedom. A chi-squared random variable 
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becomes relatively less variable as its degrees of freedom increase, 
As before, because of the required normalizations, it can be shown 
that KS CF)/S, 4) is distributed as chi-squared with k degrees of 


freedom. k is calculated (78) from 


A 


np 21 EMS esCh)e|- 


var (8, (£)} (3.41) 


Deets worth no tine that. Gt Sx FD 1S .SUbsoi buted, for Sx) in equation 
& ; = Due 2 as eae 
3.41 then k=2 since [E 18, (£1 ] Vaz {S_- Es So ED This: is 


consistent with the previous discussion. 


“aw 
“aw 


Knowledge of the approximate sampling distribution for Sf) 
is invaluable in that it enables one to make probabilistic statements 
concerning the accuracy of the estimates. That is, it is possible to 
construct intervals of confidence around the estimate SG) which will 
enclose the true value S ee (FD on 100(1-0)% of occasions on average. As 
a consequence of the chi-squared distribution of ks (£)/8_ () it can be 
shown (79) that the 100(1-0)% confidence interval for S x #) is given by 


Kone) KS (£) 
ee ao (3.42) 
z, {1-5} Ze {5} 


where 2 i and 2, Bias are the a/2 and l-a/2 points on the cumula- 


tive chi-squared distribution with k degrees of freedom. 
3.4 Consistent Cross Spectrum Estimation 
The sample cross spectrum is given by 
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R (t) is the sample cross covariance function and X, CE) and Yp Cf) are 
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the Fourier transforms of x(t) and yp Ct) respectively. Unlike ys 
which is always real, oo is usually a complex function of frequency. 

The sample cross spectrum is statistically similar to the 
sample power spectrum. The values of Sines. do not tend to a limiting 
value in any statistical sense as the record length T tends to infinity. 
Hence, See. like eae is an inconsistent estimator. Smoothing 
techniques are therefore again required. 

| The indirect approach to consistent cross spectral estimation 

involves calculating the cross covariance function between the processes 
out to lag M, where M<<T. The result is quadratically modified and 
Fourter transformed. The direct approach requires that the data records 
be segmented into subrecords of length L, with L<<T. Each subrecord is 
Fourier transformed and the products x, *(f)¥, (£) formed. 8 ny (E) is 
then obtained by eugtseLRE the contribution of cre tene from all of 
the subrecords. The eee obtained by either of these routes are 
consistent estimates of Sy SE) 

The expected value of the smoothed spectral estimator can be 


shown to be (79) 
: Dales ° f-2) d 3.44 
B ts (£)} = J, We)S,,(f-8) de (3.44) 
The bias in Say CE) will be small if eS is slowly varying over the 


spectral window W(f). 


In complex notation Sen? can be written as 


Se) = & (£) elfxy!? 6:45) 
xy xy 
where 
: = |§ (3.46) 
us (=. | ee (£) | 
and 
a -l1 IMS £ 
ec) can ) (3.47) 
al RE ae £) 


with IM and RE representing the imaginary and real parts of their 
arguments respectively. If L=2M and there is no linear or quadratic 


modification then it can be shown (82) that 


var {H, (6) ~ a chee) Sat (£) Coreates) (3.48) 
and 
TPP omega dnl ai) (3.49) 
ee 2T ye) 


These variances can of course be reduced by, for example, using the 
Hanning lag window or overlapping the subrecords depending on which 
ae. is used. It is worth while noting that the variability of 

— and Ps ey is a function of the true squared coherency between 
ad two pene Therefore, although the variances can be controlled 
to some extent by smoothing they may be dominated by the uncontrollable 
influence of the coherency spectrum. | 

A more useful function than the magnitude of the cross spectrum 
(f) is the gain function Pade Avconsistent estimator forethis 


function can be obtained from 
pg ea) (3.50) 


It can be shown (79) that approximate confidence intervals for 


a8 and Ei are given by 
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f., 1m {l-a0} is the l-a point on the cumulative Fisher's F distribution 
> a 
with 2 and k-2 degrees of freedom. ene is the smoothed coherency 


Spectrum estimator. 


iSeuD Consistent Coherency Spectrum Estimation 


The sample squared coherency spectrum is defined as 
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Mra ae aS 
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regardless of the relationship between the underlying processes. The 
sample coherency spectrum is therefore useless as an estimator of the 
true coherency between two processes. The smoothed squared coherency 


estimator 


is, (6)? 
ee (3.54) 


Ab z 
Neg) cee 
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is a considerably better estimator of Yeey E> the true squared coherency, 
and is the estimator most commonly used for this purpose. 

Using Monte Carlo methods, Foster and Guinzy (86) and later 
Benignus (87-88) have shown that ee is actually a biased estimator of 
y2_. This bias increases as the true coherency tends toward zero and 


can be considerable if the number of degrees of freedom associated with 


R 


the estimates § pe See) and = (£) are low. Using the distribution 
x 


function developed by Goodman (89) for 3 between two Gaussian 
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processes, Foster and Guinzy have shown that their maximum likelihood 

estimator ior Vey produces consistently lower estimates of y?2 than 
XY 

& 

y2 

< 


freedom associated with the estimates is less than 50. Benignus (87-88), 


a This is especially true at low im when the number of degrees of 


on the other hand, has published a set of curves which describe the bias 
in ©) As dai-hunction) of Very F) and k, the number of degrees of 
freedom. These curves show that as Vey (f) approaches zero the bias in 
v2 (f) increases monotonically. For k=4,8,16,32 and 64 the bias in 
ee ©) Boproaches: ln llnt e055. sO2ocand (01 Mesnectively.g eal tebismaliso 
shown that these curves are not sensitive to the amplitude distributions 
of the underlying processes. Jenkins and Watts (79) have shown that 
the bias in Mee) depends largely on the sum of two terms. The first 
of these becomes small as the number of degrees of freedom associated 
withthe estimates increases. The second term contains the true phase 
function between the processes. If there is a considerable delay or 
phase lag between the two processes this second term can contribute 
significant bias to ¥zy (E). 

The variance of the smoothed squared coherency estimator can 


be shown to be (79) 


Vamp neues Sahay oye Taye) 2 (3.55) 
and (82) 
ae (yay Cap ee ott aye) (3.56) 


An approximate confidence interval for sean can be constructed 


by assuming that (76,79) 
Gace) Se (3.57) 


is Normally distributed. Hence, confidence TRcervetS LOL Healey are 


given by 


=33- 


40% 


sa tay sc at 


‘i 


a 


a roe ‘art ye 


U (£) +m {l-a} (3.58) 


Ovi. a ; : P ; ; ; 
where m 1s} is “the -l= point on the cumulative unit normal distribution. 


In section 3.1, the invariance property of Vey (FD under linear 
filtering operations was discussed. Foster and Guinzy (86) have Shean that 
ve) does not necessarily share this property. If u(t) and v(t) 
are derived from x(t) and y(t) by passing the latter processes through 
the same linear filter then the invariance property extends to the 
smoothed estimators. If they are not passed through the same filter 
the invariance property only holds if the characteristics of the 
respective filters are slowly varying over the equivalent spectral 
window. If this is not the case, Hen) will tend to be on the low side 


of Vey Cf) 


3.6 The Fast Fourier Transform 


Today, the most practical route by far to spectral estimation 
is with a digital computer. Continuous data is therefore replaced by 
discrete or sampled data and the continuous Fourier transform gives way 
to the discrete Fourier transform. The discrete Fourier transform can 


be written as 


f 1 N-l -—j2Tig/N 
GS) Ba qZ0 x(q) e J q/ 


(3.59) 
where x(q) is the discrete time series and X(i) is the discrete frequency 
series. N is the number of data points, j=v-1, L=O515 2-0-2. in general 
both x(q) and XG) are complex series. To calculate all N Fourier 
coefficients the obvious way, approximately N2 arithmetic operations 

are required. In 1965 Cooley and Tukey (26) showed that if equation 

3.59 was judiciously factored, the number of operations could be reduced 
to 2NlogoN. However, N must be chosen to be an integral power of 2. For 


N=1024 this represents a computational reduction of more than 200 to 1. 
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A more striking example is given by Cochran et al. (90) who report that 
when conventional procedures are used with N=8192, the computations on an 
IBM 7094 computer require about half~-an-hour. With the Fast Fourier 
Transform about five seconds are required. On the IBM 360/67 used for the 
computations in this thesis, the calculation of all 8192 Fourier 
coefficients is reported to require about 2.9 seconds, (91).- Many 
descriptions of the Cooley-Tukey algorithm are now in the literature 
(85,90). 

When the Fast Fourier Transform is applied directly to obtain 
estimates of the power spectrum of a process, the autocovariance function 
is not a byproduct of the operation as it is with the Blackman-Tukey 
method. Jenkins and Watts(79) have pointed out that the covariance 
function can sometimes be used to advantage in reducing the bias 
especially in cross spectral estimates. Bingham et al. (80) have shown 
that to obtain the first 500 terms of the discrete covariance function 
based on 4096-500 data points, it is computationally 12 times faster to 
first go the direct route to the power spectrum with the Fast Fourier 
Transform. This frequency series is then Fast Fourier re-Transformed 
to arrive at the covariance. The 3596 original data points are padded 
with 500°zeros*to avoid circular convolutions. Before the Fast Fourier 
Transform, the most efficient route to the power spectrum was via the 
covariance function; since the Fast Fourier Transform, the most efficient 


route to the covariance function is via the power spectrum. 


—-35- 


t 
+ 


al 


i : / i’ - : a A, : ; 
eerie S 
: hr aay ; iy ' 
7 : ‘ ¥ 
7 on) Gabe nT i Ay 
Pak a 
. i 
~ lane fin / ‘ y 
var. 7 SS 2 0 
9h te he 


ARN LGR WAT ey noi weet oD _Fieaititinanatts 


‘ 


a 
, Ps. 74 Deere ee, ro 2h iste Sarr oe a 
ae c i : ° ’ 
be AN ne 7 


. ie! RS # De 
« 
a 
n é 4 e i ‘ Py rons ies 
! 
¥ » 
ie ¥ 
lan 
— - * > e 
a a 
Fh . ay * 
4 
> / j a 
g an r t 
j } j a” : 
{ . P ’ 
> = 
he an : 
i mn 
a ; ie i ) 
. “ bist { 
Se twepe oes . q 
« 7 Jf F s + Pr - 
. ? 3 eo eL 5 ead a : ipa 
we 

j ~ = é a = TG ’ 
i 


us Pinu wah, 

Tani : 
1 U2.) idee, ae ey mt, ao ; 
} Wie 


et ob he it ial ets ee LA if era pe 


CHAPTER 4 
PROCEDURE 


Bol The Preparation 


In all experiments, muscle spindles from the extensor digitorum 
longus IV muscle of the male toad Xenopus laevis were used. Often more 
simply referred to as the long toe extensor, this muscle is located in 
the hind limbs of the animal and is supplied by branches of nervus 
peronus lateralis. In the adult, the relaxed in situ length of the muscle 
(L,) is between 19 and 25 mm. The muscle contains about 50 muscle fibers 
with diameters varying from 5 to 100u. The muscle was removed from the 
limb with about 1 1/2 inches of the innervating nerve intact and placed 
in a bath containing Ringer's solution. The composition of the Ringer 
has been given elsewhere (48). The capsular portion of a muscle spindle 
was located with a microscope by shining a light through the isolated 
muscle. This region was then denuded by cutting and removing all other 
fibers from the area. Ail nerve fibers, except the single afferent axon 
which innervated the denuded capsule, were cut. The isolated preparation 
is illustrated diagrammatically near the top of figure 4.1. In the bath, 
the muscle was held taut by two wires each hooked into one of its tendons, 
Sensory discharges were obtained by measuring the potential changes 
produced by longitudinal action currents across a petroleum jelly seal 
placed around the nerve trunk. The time course of the action current is 
approximately proportional to the rate of change of the nerve action 
potential. Recordings were obtained from a total of 17 apparently 
normally functioning receptors. The temperature of the bathing Ringer's 
solution was maintained at 20°C + .5°C throughout by continuous perfusion 


from a temperature controlled reservoir. 


G22 DataiGollection 


The afferent discharges of the spindle were recorded in response 
to two classes of stimuli. The first of these could be called static 


stimuli and consisted of maintaining the muscle at various constant 
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THE EXPERIMENTAL AND DATA COLLECTION PROCEDURES. 
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lengths ranging from well below Lo Pouce the point ofulicacture, The 
other class of stimuli was of the dynamic type; upon each static stimulus 
was superimposed a small and usually random perturbation the course of 
which was recorded simultaneously with the sensory discharges. 

The wire rods hooked to each muscle tendon were securely 
fastened to the shafts of a differential puller. The puller was 
designed so that the displacements of each of the shafts accurately 
followed the voltage level at the input of the control unit. The control 
unit contained a continuously-variable signal-splitting circuit. This 
allowed adjustment of the magnitude of the displacement at each of the 
puller shafts while keeping the overall magnitude of the displacement 
between the shafts constant. It was therefore possible to keep any point 
along the length of the muscle stationary while the overall length was 
changed. The splitter control was always adjusted so that the spindle 
capsule remained stationary during the application of a stimulus. The 
length changes of the reticular zone could therefore be observed with 
a microscope while applying controlled length changes to the ends of the 
muscle. Stationarity of the capsular region also ensured that no strain 
would be placed on the nerve-capsule junction during stimulation. The 
component and circuit diagram of the differential puller along with some 
performance curves are given in appendix l. 

During all phases of the experiments, the receptor discharges 
were continuously monitored using both visual and auditory means. As 
a consequence, and through experience, a normal preparation could be 
immediately recognized. During an experiment, the audio system seemed 
to provide the first noticable sign of forthcoming abnormalities. This 
occurred as either a change in the tone of the individual discharges or 
as a general change in the pattern in which they occurred. The output 
of the preamplifier was also fed into a pulse-height discriminator- 
instantaneous frequency meter combination. The pulse~height discriminator 
is a device that is electronically capable of selecting both the lowest 
and the highest amplitude spike at its input that will be processed. It 
responds to a spike of the selected amplitude by outputting a pulse 


200 psec long and 10 volts high coincident with the fall of the input 
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spike. Responses of the spindle were recorded by first passing the 
preamplified signals through the pulse~height discriminator. This 
procedure removed the base-line noise present on the raw signal and 
resulted in clean homogeneous recordings from all experiments. The 
circuit diagram of the pulse-height discriminator is shown in appendix 
2. The instantaneous frequency meter is an analog device which computes 
the reciprocal of the time that has elapsed between successive pulses 

at the output of the pulse-height discriminator. During the interval 
between pulses, the output of the instantaneous frequency meter is 
constant and directly proportional to the reciprocal of the time interval 
between the last two pulses. The output of this device, displayed on a 
strip-chart recorder, was invaluable in, for example, determining when 
the spindle response to a change in static stimulus had reached a steady 
state. A description of the instantaneous frequency meter has been 
published elsewhere (92). 

The dynamic input-output properties of the muscle spindle 
were examined by applying random changes in length to the isolated 
muscle. The signals used as input to the control unit were obtained from 
a General Radio Company type 1390-B random noise generator. The output 
of this unit was passed through a Krohn-Hite Model 330M band-pass filter 
and recorded on an fm channel of a Precision Instrument PI-6200 tape 
recorder at 37.5 inches per second. The lower cut-off frequency of the 
band-pass filter was set at the .2 cps minimum while the upper cut off was 
set at either 200, 800 or 1500 cps. The signals recorded on the tape 
were played back at .375 inches per second resulting in random signals 
with approximately 2, 8 and 15 cps upper frequency limits. The power 
spectrum of each of the signals was estimated using the methods in 
Chapter 5. Using these results, filters were constructed to produce 
from these signals another set of signals with approximately flat power 
spectra in the range from .04 cps to the upper cutoff frequency. 

There were primarily two reasons for the choice of a random 
dynamic stimulus. The first was that the techniques of spectral 
analysis, especially cross spectrum analysis, are most reliable if the 


processes involved have relatively flat spectra. The shape of the power 
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spectrum of the sensory point discharges will be considered in Chapter 5. 
Second, a random stimulus is probably physiologically more real than 
other commonly used stimuli. However, since high frequency stimuli can 
also be unphysiological, the effect of changing the upper frequency limit 
of the stimulus had to be investigated. Problems which resemble the 
phase-locking phenomenon sometimes present with sinusoidal stimuli, were 
not encountered in any of the experiments. The time course of the length 
changes applied to the muscle were recorded by a Precision Instrument 
PI-6200 simultaneously with the muscle spindle discharges. The stimulus 
was recorded on an fm channel while the discharges were recorded on a 


Girece channel. 


4.3 Data Manipulation 


The data on the analog tape was analyzed as illustrated in 
figure 4.2. The role of the Hewlett Packard 2116B computer was to 
digitize the analog data and convert the resulting binary records into 
a format that could be easily and quickly read into the memory of the 
IBM 360/67 computer. The '360 was in turn programmed to perform the 
required computations and to produce the instructions for the Calcomp 
PLOtEer UNLE. 

The spindle discharges recorded on the direct channel of the 
analog tape recorder were preprocessed in one of two ways prior to 
analog-to-digital conversion. This, in both cases, consisted of 
reshaping the spikes into narrow square pulses of about 200 usec in 
duration and 10 volts in amplitude. These in turn were either used to 
trigger a bistable multivibrator (flip-flop) or passed through a low- 
pass electrical filter. The waveform presented to the analog to 
digital converter was then either a square wave whose state changed with 
the occurrence of each action potential or a waveform resembling electronic 
shot noise (93). Preprocessing of the signal representing the changes 
in the muscle length was confined to a gain adjustment of the waveform 
so that maximum resolution over the +1 volt range of the analog to 


digital converter could be obtained. 
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Each analog to digital conversion resulted in the production 
of a 10 bit binary number. The +1 volt input range of the converter is 
therefore quantized into 1024 equal increments. According to Bendat and 
Piersol"(76). the effect of quantization upon the continuous waveform can 
be considered as an additive rms noise on the original waveform. The 
results of these authors can be used to show that the additive rms 
noise due to the 10 bit quantization of the +1 volt range of the converter 
amounts to about .3 millivolts rms. 

The 2116B computer was programmed to accept binary numbers 
from the converter at regular intervals until a total of 8192 numbers 
had been collected. Each number was then converted into 6-level binary 
coded decimal (BCD) characters and written on a 7-track digital tape. 
All translations from BCD to the internal code of the '360 were performed 
by the '360 during the execution of the ‘READ* command in the Fortran 


IV source language. 
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CHAPTER 5 
ANALYTICAL METHODS 


Bela Gataks CLCaL Techniques 


If the flip-flop output of the pulse shaper is sampled 
regularly, the occurrence of an action potential can be detected by 
testing for the condition where two consecutive samples are not equal. 
This procedure was used to obtain the information needed to construct 
pulse-interval histograms of the spindle responses to static stimulation. 
The instantaneous frequency meter indicated that under all conditions of 
static stimulation, the discharge rate of a normal preparation seldom 
if ever exceeded 50 pulses per second (pps) and always averaged less 
than 20 pps. - The sampling rate for the analysis of interval statistics 
was chosen to be 200 samples per second. The resultant bin width for 
the interval histograms was therefore 5 msec which divided the average 
interpulse interval into at least 10 units and resulted in histograms 
of useful stability. 

The histogram of intervals between successive events (histogram 
of first-order intervals) was used to estimate the two parameters of the 


Gamma density 
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where g and 8 are the required parameters. [I(g) is the gamma function 
of g. g is often called the order of the Gamma density. p(t) is the 
probability density of an interval t. The mean and variance of this 


density are 


and 
E{(t-t)2} = var{t} = 0,2 = 8, (5.3) 


: its. ? 

{ ar #s A ree. 
win? te oe ae Pie” “a 4 Rosi yes' oe 3 Bs ay stay etal ot ‘omy he ne te 
=. ae 

5 9 4 

LO OES Lapin se IE Bieber Siu pert Say a gas feeds: ERO ET ‘eran “s 

. . 7 r F 
Web Tes (CR we iia Teen a A re cee oa ee aiid nat 


P . ‘ae be , ’ ae 7 ra Saat 
erie Aa ede tai th etaa lac (oni) Teta eae a eae tea | 


244 10) Bee Peuies Teer. east). ae a Fil} Lae oO: bniaas AAJA 0 (etawaaicd | 


ade Aw Cee ee 
1 i i " \ W 
oe 
sa mit eee As 
al iy She ei" 
Af dt Ss fe ante. 
a » r Fa ] Ph ie 7 
“h 
7 hex = ca e - a] eS 
DOES. C2 eet ee ey RA Si, tired ait: cv rf sit at : 
Veo Hah, ed ae EES, rte ©): Pa ec eerie oc? Ue a ba cody I ahp tt 
—* ba 
{ Hs unre eels Bale <2 rc eRe Thi 1 OGL ast! bhi ; nara ln 
. oa ; Ls 


POTS HON OS? Faas a ae een el end | RY ebran Ti 1 nite 4 


fit Shyness Gy tedh (agen Dit Saye jie ite} Gx iw 
. a") 

OT oar wid SRR te) See |. Caen 1 eu ae Woe. i whe : 

SR IAT) 202 Hen eke oie) pe sare ir eas Ta: Raa gles ee 

BARTS eT aed PGS) ‘lah. sOuiagy Ul Peete bce fps 

i; ; | oe bhi 8 

Oo etilh eeppue HBR oid Day ed ai Lewy = i ia ‘Glin’ bit ae 


in : , ‘ald ; it ae ix i oF ). . 
i" ) i a sete { : Dy a ie n 
a ‘ cals ES 7 : a Se i} #6 nsf ris) ' 

Y \ 


ie iN Aiton “y's Wis 
bo PA yby® ‘cits soph) Cae expctialiek ig Bi pad tp ae ans 5 | 


ould ei C4) 9 ot Hest sas ate te tht pt bait tai mbes 
etd 4) howe layiy, sit ee be er, ‘ae a es 


respectively. 

Estimates of 8 and g can therefore be obtained by using the 
mean and variance of the intervals. However, Cox and Lewis (94) have 
shown that only poor estimates of these parameters are obtainable by 
using the sample mean and variance. This is especially true for low 
values of g (g<10). Estimates of 8 and g using the sample mean and 
variance are called inefficient estimates. Bendat and Piersol (76) 
define an efficient estimator as that which produces the least mean 
Square error between the parameter itself and its estimator. 

The class of maximum likelihood estimates is known to exhibit 
the minimum variance that can be achieved by any unbiased estimator (79). 
For large data samples this estimator is approximately unbiased and 
Normally distributed. The maximum likelihood estimates for 8 and g 
are obtained by maximizing the likelihood function for the N observed 
times between the events with respect to these parameters. The likelihood 
function is 


N 
Ne ,N , g-l, -8) t, 
B Ca ie )e fay - 


LCG, 256) = 1 


[r(g)]” 


where t. refers to the ith interval. The results of such an analysis 
a 
(94) indicate that the maximum likelihood estimator for g is the solution 


of 


‘ N 
log (g) - ¥(g) = log t - } log t, (5.5) 
e e : ei 
i=l 
N 
where t is the sample mean given by 
a N 
So (5.6) 
i=l 
N 


~ he 


Bs a hare a any Lik acl ha Ok atin teh C 
ih ! ra 
: 
r ve he ) 
yy } ’ a I My i" : 
4 4 a ae / ' 5 i 
J { cy ‘ ' 
Al’ j 7 | ‘ . 
i ' , A , ‘ ] 4 2 
lee a ; 
(rota }) | 4 + t 5 as o> 
F ’ 
. 3 Si) OS 
‘ 
‘ ny.) *) 
- r z .$ - 4 
i 7 ES 
i | ji : 
" 2 
hs : 
x ra . 
. Oe th. 5 j 
‘ f 
‘ ia 
f f a 
' P , 4 ‘ 
: ; UPN ES 
' 5 a ; 7 \ 
- tr ee a, esti i 
; : 
6 : n eal 
' on} ; | 4 
4 4 erik \ ieee as 
g \ ; 
a 4 a “Sa a 
av] ¢ fi > i 
| ‘ae 
fl : y 
al ; i , i, Be ” 
‘ - a we. { 
SEReAens. HR ge oo eae et, BE “ “2 — atitaat an Sed mai: mi eR 
rectal tae Pe. o ok 3q7 Git Sc 3B val bolt a qt i0o kee 8 
\, oe) ) r r a 
oh v1 ; ’ : ee ee i i 
ie one | i a en : 
{ \ ry 7 of Cee | 
W ace nent 
i 7 a ¥ A" o- ‘i F are irae 
c) M . yh re bP 2 bi ~~ @ J ma ee $e. 
ant seg 0 fy), AU Peele we 9 is oF ii, 7 a ane 
rr ee ee a ers Ube. 2S | : hy ssl OR aeer th 
: oe | 7) Tide an , A Aue bal te Ur “te oe 7 ve a hee nh . 
ee et a an ee 
7 a t= Pr! Lat : i eo 7h i 
4 Sie] Bi see ot =f ie mia" Mya 7 eae: 
i fae oS cal ne ia Ny 
iy, i i) ee, AY 7 » T vin ca i 
f ; sia >. ie i se 
> 4 Ae i Ad rw 7 ay iat 


: ys oe, : im eae a 1 ' 
ip ee: iy 7 ae: 
i =e) ‘i re 

wre 


eae o Leahy 


ua , lh ‘, ey 
2 “ Y 


¥(g) is the so-called digamma function and is defined as 


2 a, 
Y(g) = d log f(g) 


se . (ana) 
dg 
a N 
The quantities t and : log t, are therefore sufficient statistics for 
i=1 


maximum likelihood estimates of B and g. Using these, equation 5.5 

was solved iteratively. The initial estimate for g was obtained from 
the sample mean and variance. To evaluate equation 5.7 subroutine 
DLGAM which calculates log T(z) to 16 significant figures was used. The 
digamma function was obtained by replacing the differential operator 


d/dg in equation 5./7 by the finite difference operator 6 where 
Sense ti Gh erana regmliteaed Agy leet (ato =) (5.8) 
"e cee 2 e Zz 


It can be shown (95) that the relationship between the differential 
operator and the finite difference operator is given by 


a 


= <2 sinh”? (5.9) 


_d 
dg 


where Ag is an increment in g. The first few terms of a Taylor's series 
expansion of the right hand side of equation 5.9 were used for equation 


5.7. The maximum likelihood estimate for 8 was then obtained from 


(57240) 
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The Gamma density with the maximum likelihood parameters was scaled to 
the same area as the first-order interval histogram from which it was 


estimated and drawn on the same axis. 


The characteristic function of a random variable t is defined 


as 
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ge ty and ty are two independent random variables with a common 
distribution function then the characteristic function of their sum is by 


definition (94) 
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The characteristic function of the sum of three or more independent 
random variables with common distributions follows from equation 5.12. 


The characteristic function of a Gamma density can be shown to be (94) 
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‘ 2nfy8 | (5.13) 
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The sum of two independent random variables each distributed with 
identical Gamma densities is therefore also distributed with a Gamma 
density. It follows from equation 5.12 that the parameters of the sum 
distribution are 8 and 2g. Parameters of the sum of three random 
variables would be 8 and 3g and so on. 

As a test for the independence of successive interpulse inter- 
vals, the parameters of the maximum likelihood Gamma function of the first- 
order spike intervals were modified as required to obtain the density 
functions which describe the sum of two and three successive intervals 
(second and third-order intervals) assuming all intervals are independent. 
These curves, properly scaled as to area, were drawn on the same axis 
as the second and third-order interval histograms obtained from the 
data. 


It can be shown (96), that the power spectrum of a point 


process whose intervals exhibit the properties of a random variable, is 
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where ees) is the characteristic function of the ith order interval. 


abe 


A is the area or strength of the points in the process, If successive 


intervals are independent then as a consequence of equation 5.12 (see 


(96)) 
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eT 
s,@) = ice (5.15) 
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Re represents the real part of the complex quantity in the brackets 
which follow it. The power spectrum of an independent Gamma process is 


then 


s (ey = ACB fitore 1 ] (5.16) 
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B 


This equation was programmed on the computer and values of eee, were 
calculated at various frequencies for a number of different g. To test 
the sensitivity of Eeoue to small changes in the shape and the form 

of the underlying density function, the power spectrum of a high order 
Gamma density was compared to that of a Gaussian density both with the 
same peak and variance. The characteristic function of a Gaussian 


density is 


where t and 2 are the mean and variance of the intervals respectively. 
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The power spectrum of the independent point process with this distribution 
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The serial correlation coefficient between successive intervals 


is defined as 


(t,~t) (t,.,-t)} 


Sd ge ay ee em ee 
ire a (5.19) 
where Pn? h=0,1,2,... is the serial correlation between intervals 


separated by h events. If the intervals are drawn independently from a 
common distribution then the expected value of all the coefficients is 
zero. Since serial correlation coefficients are one of the most powerful 
tests for the serial independence of intervals and also one of the easiest 
to implement, the results of such an analysis on the spindle discharges 
have been included. Estimation of Ph requires that t and oe be estimated 
from the sample data. To avoid the bias introduced by the use of the 
sample mean and variance, the formulas suggested by Cox and Lewis (94) 
were used. These formulas are repeated here. The estimator of PL Ls 


given by 


h (5.20) 
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A ] 
For N large and p79; pata) 4 has an approximately unit Normal 
distribution. Therefore, the hypothesis of independence can be rejected 


at the a level if 


lon! Chery | (5.24) 


where Bits ae) is the l-a/2 point on the cumulative normal distribution. 
To reduce the possibilities of obtaining false positive serial correla- 
tions that result from long-term trends or nonstationarities in the 
process (2), the data was always divided into a number of segments. The 
first twenty serial correlation coefficients were calculated using 
equations 5.20 to 5.23 for each segment. The final estimate for the 
twenty coefficients was obtained by averaging the contributions from 


each segment. 


5.2 Sampling the Point Process 


If the spindle discharges are treated as a series of point 
events (Dirac delta functions) the sample Fourier transform of the data 


is given by 
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where 6 (t-t,) is a Dirac!delta function at t=t.. When treated this way 
the estimation of the power spectrum of nerve spas is extremely well 
suited to digital computer evaluation since only the summation of cosine 
and sine terms at the instants in time at which the spikes occur is 
required. Knox (97) has used this method to obtain estimates of the 
power spectrum of the discharges from receptors in the statocyst organ of 
the lobster. A serious shortcoming of this method is apparent from the 
behavior of K(f) at £=0. Xn Cf) at this point is a spike of amplitude 

N, the number of action potentials in the sample data. Since X,,(0) 
stands far above x, (CE) for £40, leakage through the spectral window will 
introduce considerable bias into the estimates of the power spectrum near 
zero frequency. 

If the method used to detect the occurrence of events is by 
regularly sampling the flip-flop output of the pulse shaper in figure 
4.2, it is necessary that the sampling rate be much higher than the maxi- 
mum frequency component of interest. The errors resulting in a sampling 
rate that is too low would probably be somewhat similar to the quanti- 
zation error which results when a continuous waveform is sampled with too 
few bits available for the conversion. 

It is a well known result, that if a continuous waveform is 
sampled regularily at a rate which is at least twice as fast as the 
highest frequency component in that waveform, then the sampled version of 
the continuous waveform contains an undistorted copy of the power spectrum 
of the original waveform (75,98). If the signal contains components 
higher than half the sampling frequency a so-called aliasing of frequencies 
occurs. In this case, the signal must be sampled at a higher rate or, 
the undesirable components removed before sampling by filtering. Half 
the sampling frequency is called the Nyquist or folding frequency. 


Aliasing of frequencies refers to the phenomenon whereby power in a 
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waveform above the folding frequency masquerades as power at frequencies 
below the folding frequency. The general rule for the aliasing of 
frequencies has been given by Bendat and Piersol (76). For any frequency 
f in the range OSfSE where fo is the Nyquist frequency, the higher 


frequencies which are aliased with f are defined by the infinite series 
(2f +f), (4f +f), (6f +f), ae (3.26) 


Once a signal is sampled, nothing short of resampling after either 
filtering the signal or increasing the sampling rate can possibly remove 
any aliasing which may have occurred. 

The estimates of the various spectra which involved the point 
discharges of the spindle were obtained, with one exception, by filtering 
the shaped discharges. Equation 5.25 was applied in one instance by 
sampling the flip-flop output of the pulse shaper at a rate much higher 
than the highest frequency of interest. In all other cases, the shaped 
pulses were passed through a low-pass filter. The filter was designed 
so that after estimating the spectrum of the filtered pulses, the effect 
of the filter could be removed from the estimates up to 80% of the folding 
frequency with negligible error due to aliasing. Filtering the pulses 
has the effect of converting the point process into a continuous waveform 
so that all the ideas of chapter 3 are directly applicable. Nelsen (96) 
has shown that if the interevent intervals of a point process are drawn 
independently from a common density function, then, if this density function 
is bounded, the power spectrum of the process tends to a constant value 
with increasing frequency. There is evidence that interspike intervals 
from amphibian muscle spindles are independent (54,72) and that the 
underlying density function is bounded (54,71). Work with equations 
5.16 and 5.18 (see chapter 5) which are based on these assumptions, has 
shown that the power spectrum of a point process with about as much 
variability as the discharges from an amphibian muscle spindle is flat in 
the reeiion of the spectrum above the mean discharge frequency. The sampling 
rate foe spectral estimates was chosen at 100 samples per second, The 


power spectrum of the sensory discharges was assumed to be flat in the 
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region above the folding frequency of 50 cps since the average discharge 
rate was always well below this. Also, no stimulus was applied to the 
spindle which contained any significant power above 30 cps. Therefore, 
because of the filter, reliable estimates of the various spectra could 
be obtained to about 40 cps. The band of frequencies topped by 40 cps 
seemed a physiologically relevant choice. 

Consideration of the bias in the estimates of the various 
Spectra seemed to indicate that the filter characteristic in the passband 
should be reasonably flat and should not introduce excessive phase lags 
between the input to, and the measured output from, the spindle. The 


filter used had the following magnitude and phase characteristics 


respectively, 
MGese ee (5.27) 
ae 
1+i— 
\2 
and 
r Y2£ > 
B= —tan 9 wor (5.28) 
Leper 
ZL 


The circuit diagram of ‘this filter is given iniappendix 8., Li, the. power 
spectrum of the unfiltered spikes is Sop then the power present in 


the signal after filtering will be given by 


Sr Lhe ee Ws aa Ga ie (5.29) 


using equation 3.7. If the power in oom at 40 and 60 cps was approxi- 
mately the same, then the power in SF) at 60 cps would be down to 

less than 2.5% of that at 40 cps. At 70 cps it would be down to less 
than .2% of that at 30 cps. Since 50 cps is the folding frequency, the 


aliased power in S__(f) at 40 cps from 60 cps and higher frequencies 
ax 
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would be about 2.5%. At 30 cps aliased power would be less than .2%. 
Therefore, assuming that ace is reasonably constant above 50 cps, 

and this can be established by plotting SE)» estimates of Spo 

can be obtained up to at least 40 cps with little error due to aliasing. 


Estimates of SPBoe from 53) are obtained using equation 5.29. 


3.3 Methods of Spectral Estimation 


Spectral estimates were obtained using the direct approach and 
the Fast Fourier Transform subroutine PS301A (91). The method was that 
suggested by Welch (83), slightly modified. The data, both single and 
double channel (spindle under static and dynamic stimulus respectively), 
was sampled 100 times per second. Exactly 163.84 seconds of data were 
used for each analysis. Since the HP2116B, as programmed, had capacity 
for 8192 data points, power spectrum data was stored on the digital 
tape in two blocks of 81.92 seconds each. Data for cross-spectral 
estimates was stored in four blocks of 40.96 seconds each. Welch (83) 
suggests that in order to obtain a near maximum reduction in the variance 
of spectral estimates from a fixed amount of data, the data segments from 
which the sample spectra are obtained should be overlapped to about one- 
half their length. This suggestion was followed. The data was divided 
into segments of 1024 data points. Sample spectra were obtained from 
POivts i. bo L024, points o12, to VasG, LOZ24ebogz04s andiso om Untilval L 
points in the data block were exhausted. The procedure was repeated for 
all blocks and the spectral estimator obtained by averaging the contri- 
butions from all the sample spectra (See equations 3.34, 3.35 and 3.43). 
The power-spectrum estimator was therefore the result of averaging 30 
sample power spectra while the cross-spectrum estimator was the result of 
averaging 28 sample cross-spectra. 

Welch (83) has shown that the variance of the power-spectrum 


estimator described in the preceding paragraph is 


142 } oe E(i) (5.30) 
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where 


rt A 
| : w(h) vont) 
: L hal 
OD = i 5 (5.31) 
la 
h=1 


K is the number of data segments in the estimate, L is the number of 
points in each segment and D is the number of points between the starting 
points of the segments. For the power spectrum estimator used in this 
thesis K=30, L=1024 and D=512. wh) is the data window. It is of interest 
to note that if D2L, €(i)=0 for all i and hence var(S_(E)} does not 
depend upon the shape of the data window. Sloane (81) has shown however, 
that windowing the data does result in a loss of spectral resolution. 

If the data segments are overlapped, the shape of the data 
window does effect the variance of the estimator. It can be shown, 
using equations 5.30 and 5.31, that overlapping the data segments by 
half their length with a rectangular data window results in a variance 
reduction to 75% of that if the segments are not overlapped. If the 
Hanning data window is used then 50% overlap of the data segments results 
in a reduction in variance to about 55% of that when the segments are 
not overlapped and the rectangular data window is used. All spectral 
estimates in this work were obtained using the Hanning data window. There- 
fore, using equations 3.41, 5.30 and 5.31 it can be shown that the power- 
spectrum estimator had 56 degrees of freedom while the cross-spectrum 
estimator had 52 degrees of freedom. Confidence bands for estimates of 
the power spectrum, gain and phase functions and the squared coherency 
spectrum were thus obtained using equations 3.42, 3.51, 3.52 and 3.58 
respectively. The statistical tables given hy Bendat and Piersol (76) 
were used. Gain and phase estimators were corrected for the presence of 
the spike filter by using equation 3.6. 

It can be shown (81) that the spectral window resulting from 


the Hanning data window is given by 
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This. equation is shown plotted jn figure 5.1 with a normalized amplitude 
of one. Leakage problems were all but eliminated from the estimates by 


this spectral window. 


5.4 Plotting the Results 


Thep Fast (Fourier Transterm of 1024, dataspoints results, dn 
1024 Fourier coefficients. These coefficients sample the frequency 
Space from zero frequency up to the sampling frequency in steps of 
1/10.24 cps. All estimators were plotted on a logarithmic frequency 
scale. This transformation of the frequency scale resulted in points that 
were widely separated at the low frequency end and highly clustered 
toward the high frequency end. To improve the low frequency resolution 
all plots shown in the results were obtained by padding the data segments 
out to 2048 points with 1024 zeros before Fourier transformation. This 
procedure, suggested by Bergland (85), increased spectral resolution to 
1/20.48 cps or to about .05 cps and improved considerably the appearance 
of the estimators at low frequencies. Power spectrum and gain function 
estimates were plotted in decibels (db) while phase function estimates 
were plotted in degrees. 

Of the 2048 points which resulted from each analysis for each 
estimator, only the first 800 were required to span the region of interest 
from .05 toi40°¢ps. The plotting function as written’ Yequired three 
arguments: the abscissa dimension, the ordinate dimension and the 
minimum distance between the plotted points. A fourth argument required 
for cross-spectrum analysis was the maximum frequency of interest. This 
was necessary since, for example, coherency studies are meaningless at 
frequencies outside the bandwidth of the applied input. Before plotting 
points, all data was scaled to the abscissa and ordinate dimensions 
supplied. Points were plotted starting from the low frequency end. 


No point was plotted which fell closer to the previously plotted point 
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Figure: D.L 


SPECTRAL EQUIVALENT OF THE HANNING DATA WINDOW. 
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than the minimum distance parameter. This point was instead averaged 
with succeeding points until the logarithmic mean of the abscissae of 
the first and last point averaged exceeded the minimum distance para- 
meter. All 800 points were used in this way. Because of the logarithmic 
frequency scale, averaging became more and more prevalent as the plotting 
progressed. The result of this was that the statistical confidence in 
the plotted points was increased as more and more averaging was required. 
‘ne width of the spectral window shown in figure 5.1 is not 
changed by padding the data segments with zeros (85). Therefore averaging 
two adjacent points does not halve the variance of the new estimator 
since the average was formed using nonindependent constituents. It is 
not until 5 adjacent points are averaged that two independent points are 
combined. Therefore, when plotting confidence bands around the various 
estimators, the number of degrees of freedom associated with a particular 
plotted point was increased by the base value (56 for power spectra, 
52 for cross spectra) each time the number of averaged points reached 
a multiple of 5. 
Appendix 4 contains a copy of the Fortran IV source program 
used to perform the cross-spectral analysis. The program instructed 
the computer to read two channel data from the 7 track tape, to calculate 
estimates of the squared coherency spectrum and the gain and phase 
functions together with their respective confidence bands, and to plot the 


results. 
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CHAPTER 6 
RESULTS 
6.1 Static Behavior 


The results reported in this chapter were obtained from three 
separate muscle spindles. The responses of these however, were represent- 
ative of all the normally functioning preparations examined. The instan- 
taneous frequency of the discharges from a single muscle spindle under 
conditions of static stimulus is shown in figure 6.1A. Recordings of 
the discharges from this particular spindle were obtained at various 
steady lengths which ranged from about 1.5 mm less than the resting length 
of the muscle in vivo (Ly) to the point of fracture about 1.3 mm greater 
than Ly: Ly for this muscle was 23 mm. After changing the static 
stimulus, data was recorded only when the instantaneous frequency of the 
discharges had apparently reached a steady state. The occurrence of 
an action potential was detected by sampling the flip-flop output of the 
pulse shaper 200 times per second. Some results from the subsequent 
analysis are shown in table 6.1. The sample means calculated in table 
6.1 were not biased by the spike detection technique. On the other hand, 
the sample variances were biased positively by about 8 msec*. This 
bias has been removed from the sample standard deviations and the sample 


coefficients of variation given in table 6.1. 
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Muscle NOweoOL Mean Interval Standard Coef. of 


Length Intervals Interval Variance Deviation Variation 
L -1.5mm Ee ie) 66 msec 420 msec* 20.3 msec geome 
L, 2529 65 a7 oad 729 
L +0. 21am 2326 70 538 73.0 ess. 
L +0. 4mm 2370 69 443 20.8 30) 
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Table 6.1. Some Statistics of the Static, Response 
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Figure 6.1 


INSTANTANEOUS FREQUENCY OF THE DISCHARGES FROM A DE-EFFERENTED AMPHIBIAN 


MUSCLE SPINDLE IN RESPONSE TO STATIC AND RANDOM DYNAMIC LENGTH STIMULI. 
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Table 6.1 shows that the spindle responded to a static stimulus only 
slightly. It is not until the muscle was stretched to near the point of 
fracture that any significant changes are noticable in the response. The 
variability in the intervals, defined as the ratio of the rons 
deviation to the mean (coefficient of variation),was largely indifferent 
to the static stimulus. These results are in agreement with those of 
Buller et al. (54) in as much as for a 70 msec mean interval from the 
frog muscle spindle they obtained a standard deviation of 19 msec. 

First, second, and third order spike interval histograms at a 
muscle length of L +0. 2mm (table 6.1) are shown in figure 6.2. The first 
order interval histogram was fitted with the maximum likelihood Gamma 
density curve estimated from the data. The second and third order 
histograms were fitted with Gamma density curves whose parameters were 
obtained from the first order curve by assuming that all intervals were 
statistically independent. The data was fitted extremely well by the 
Gamma curves and therefore supports earlier findings related to both the 
distribution of, and the dependence between, spindle interspike intervals 
(54,71,72). The maximum likelihood estimates for the parameters of a 
number of Gamma densities obtained from the spindle intervals under 
different static stimuli are compared in table 6.2 with the estimates 
of these parameters obtained from the sample mean and variance. The 
bias was removed from the sample variance before calculating the 
parameters. 
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Table 6.2. Estimates of g and 8. 


If the premise that interspike intervals were independent 
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DISTRIBUTION OF ENTERVALS DURING STATIC LENGTH STIMULATION. 


and distributed with a Gamma density is accepted, then the power spectrum 
of the point sensory process can be calculated by using equation .6, 

The results of such a calculation are shown in figure 6.3 for Gamma 
densities of various orders, In this figure w = 2rf. All of the curves 
are normalized to the same high frequency value to enable comparison. The 
Poisson process, with its exponential distribution of intervals, 
corresponds to a Gamma density of order 1 and has a flat power spectrum 

at ald frequencies. As the order of the Gamma density is increased 

from 1, the spectrum contains less and less power at low frequencies. 

The transition between the very low frequency power and the high frequency 
power is much more abrupt and tends to occur at lower frequencies as 

the order of the density is increased. 

Estimates of the power spectrum of the spindle discharges were 
obtained initially in two ways. The first was as described in section 
5.3, that is by filtering the shaped pulses, applying a Hanning data 
window to the data segments from which the mean had been removed, and 
finally, averaging the sample spectra obtained from the modified 
segments. The results of such an analysis on a statically-stimulated 
pedi are shown in figure 6.4A. The 95% confidence band for the 
estimate is also shown. Figure 6.4B shows an estimate of the power 
spectrum of the same data obtained by sampling the flip-flop output of 
the pulse shaper 200 times per second. The continuous Fourier transform 
of each data segment was obtained using equation 5.25. Figure 6.4B is 
an excellent demonstration of distortion at the low frequency end of the 
spectrum due to leakage from the spike in the spectrum at zero frequency. 
If the distortion due to leakage is temporarily ignored, both figures 6.4A 
and 6.4B indicate a relatively flat low frequency spectrum with an 
increase of about 5db occurring in the power between 1 and 10 cps. This 
similarity is UGne ied token proof of the validity of the impulse 
filtering method of spectral estimation outlined in section 5.2. 

Figure 6.5 shows the estimate of the power spectrum of the data 
used in figure 6.2. The estimated spectrum is fairly flat in the 
region less than about 2 cps but increases by between 10 and 15 db from 


there to 10 cps. The maximum likelihood Gamma parameters for the first 
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Figure 6.3 
THEORETICAL POWER SPECTRA OF POINT PROCESSES WITH INTERVALS 


DESCRIBED BY GAMMA DENSITIES OF VARIOUS ORDERS. 
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Figurevoy5 
ESTIMATE OF THE POWER SPECTRUM OF THE POINT PROCESS 


IN FIGURE 6,2. 
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order intervals in figure 6.2 are Be 4 and 3=9.5, Figure 6.3 indicates 
that the power spectrum of a Gamma density of order 9 undergoes a change 

of slightly less than 10 db. In this respect, the theoretical and 
experimental results agree. However, figure 6.3 shows that the trans- 
ition from low-frequency power to high-frequency power for g=9 is all 

but complete at w=8. This means that the ee pecnenial spectrum should 
have reached its high frequency value by about .03 cps. This is definitely 
not indicated by figure 6.5. 

There are two possible explanations for the discrepancy between 
the theoretical and the estimated spectra. These are that the intervals 
are not distributed as a Gamma random variable and/or that the intervals 
are not independent. To test the sensitivity of the power spectrum to 
slight changes in the shape of the underlying density function, the 
theoretical spectrum of a point process whose intervals are distributed 
as a Gaussian random variable was compared with the spectrum of the 
Gamma density of figure 6.2. The Gaussian density function was chosen 
to have the same peak and variance as the Gamma density. The two density 
functions are shown in figure 6.6. Figure 6.7 shows the respective 
spectra each normalized to the same high frequency value. This figure 
shows that the power spectrum of an independent point process is not 
particularly sensitive to slight variations in the shape of the under- 
lying density function. The discrepancy between the theoretical and the 
estimated spectra is apparently not attributable to the choice of density 
function. 

Figure 6.8 shows the first nine serial correlation coefficients 
between interspike intervals from the spindle in table 6.1. The left side 
of figure 6.8 was obtained from the intervals in figure 6.2 (L +0. 2min) 
while the right side was obtained at a longer muscle length (L,t1. 1mm). 
Formulas 5.20 to 5.23’: were used with segmented data to obtain the results 
in figure 6.8. It can be shown that the spike detection technique will 
result in a slight negative serial correlation between adjacent intervals 
even if the interspike intervals are statistically independent. This 
Beane the first serial correlation coefficient amounts to about .004 


in the left side of figure 6.8 and about .007 in the right side. To 
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Figure 6.6 


THE GAUSSIAN AND GAMMA DENSITY FUNCTIONS. 
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THE THEORETICAL POWER SPECTRA OF POINT PROCESSES WITH INTERVALS 


DESCRIBED BY THE DENSITY FUNCTIONS IN FIGURE 6.6. 
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independent intervals were pe ewaead from equation 5.24. In spite of 
the bias introduced by the detection technique, it is statistically 99% 
certain that adjacent intervals were negatively correlated. That is, an 
interval of less than the mean length tended to be followed by one that 
was greater than the mean and vice versa. This tendency to oscillate 
about the mean is apparently not of sufficient strength to result in 
Significant contributions to the estimated power spectrum of the process. 
The large positive value of the fourth coefficient in the right side of 
pigure 6.0 1s Jikelysdue to. statistical fluctuations in the estimator 
Smee it 1S not present in the left side of the figure, It is to be 
expected that for purely statistical reasons, about one in a hundred of 
the coefficients will be outside the 99% confidence region. The obvious 
weakness of the serial dependence between adjacent intervals coupled with 
much shorter data samples probably accounts for the earlier results 
concerning the independence of spindle interspike intervals (54,72). 

That neighbouring intervals are correlated might be expected from what 
is known about the initiation and conduction of action potentials in 
nerve fibers. 

The power spectrum of a point process is apparently not exception- 
ally sensitive to serial correlation between the intervals. In a sense, 
use of the power spectrum to test for serial correlation between intervals 
is much like the test applied in figure 6.2. That is, both tests are 
rather qualitative in that any serial dependencies would probably 
become obvious only when a strong correlation exists. In addition, 
neither test would provide any information concerning the nature of any 
obvious correlation. But unlike the serial correlation coefficient, the 
power spectrum would probably maintain its sensitivity to serial dependence 
in the presence of trends that are long compared to the mean interval. 
These trends or nonstationarities appear as very low frequency components 
in the power spectrum. Since leakage between frequency bands during power 
spectrum estimation can be effectively controlled (see figure 5.2), the 
effects of the nonstationarities can be easily confined to. the appropriate 


portions of the spectrum. 
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6.2 Dynamic Behavior 


The effect upon the instantaneous frequency of epi ate 
discharges of three different dynamic stimuli is shown in figures 6.1B, 
C,D. Each of the stimuli is a random signal and each has a different 
upper cutoff frequency. Estimates of the power spectrum of each of the 
Signals, enclosed in 95% confidence intervals, are shown in figure 6.9. 
The signal described by figure 6.9A which shows a relatively flat 
Spectrum up to about 2cps will be referred to as input A; the signal in 
figure 6.9B with a relatively flat spectrum to about 8 cps as input B 
and figure 6.9C with a flat spectrum to about 15 cps as input C. Figures 
6.1, B,C,D which show the spindle responses to inputs A,B and C respec- 
tively, indicate that increasing the upper cutoff frequency of the 
dynamic stimulus resulted in discharges of higher instantaneous frequency. 

Figure 6.10 shows estimates of the power spectrum of the 
spindle response to the three inputs at three different lengths. The 
steady state spectral estimates are included for comparison. A1l 
confidence intervals for this and succeeding figures are shown at the 
95% level. The preparation in figure 6.10 is that described during static 
stimulation by table 6.1. The top row in figure 6.10 was obtained from 
the spindle at a static stimulus of Lj-i- 5mm, the second row at L tO. 4mm 
and the third row at Lt. Omm where Lo was about 23mm. The amplitude of 
all the inputs was about 45um rms throughout. Columns 2,3 and 4 correspond 
to inputs A,B and C respectively. 

Figure 6.10 indicates that the response of the receptor to 
the various inputs was largely an addition to the spectrum that 
existed before the dynamic stimulus was applied. This implies that the 
steady state spectrum was acting as the carrier for information concerning 
the dynamic state of the muscle length. In the band of frequencies 
studied, application of input power resulted in the addition to the 
output power of the receptor, power at only. those frequencies at which 
an appreciable input existed. Interaction between the carrier spectrum 
and the input spectrum, apparently did not result in any sideband copies 


of the input. By comparing the responses of the spindle to the same 
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dynamic stimulus at different steady-state llengths, iteis clear that the 
dynamic gain of the receptor depended strongly upon the length of the 
muscle. This is best demonstrated by input A. Here, the shape of the 
spindle response is almost identical at all three lengths with the 
exception that at the longest length the response was about 5 times 
greater (7db) than at the shortest length, The responses to inputs B 
and C are complicated by the fact that their spectra extended into the 
transition region of the carrier spectrum. Also, because their broader 
bandwidths resulted in relatively lower power densities, the power 
spectrum of the spindle response to inputs B and C was not as definitive 
as it was to input A, especially at the shorter muscle lengths. However, 
another feature of the spindle response is apparent froin inputs B and 

Cc. These inputs have approximately equal power densities in the region 
below 7 cps. The flat portion of input C extends to somewhat higher 
frequencies. In the regions of the spectrum where inputs B and C have 
nearly equal power densities, it might be expected that the power density 
of the receptor response to these inputs would also be roughly the same. 
This was definitely not the case (cf figure 6.10 row 3 columns 3 and 4) 
and significantly less power existed in the response to input C 
especially in the region 1 to 7 cps. This then is an indication that 
the spindle gain was also dependent upon the bandwidth of the length 
variations applied to the muscle. 

Estimates of the squared coherency spectrum between the spindle 
response and the impressed input are shown for the various cases of 
figure 6.10 in figure 6.11. Coherency between input and output was 
always greater with input A at all muscle lengths than with either 
input B or C. This might be expected since the power density of input 
A was so much greater in the region less than 2 cps than that of inputs 
B and C. However, inputs B and C had equal power densities to about 7 
cps and yet the spindle response tended to be more coherent in that 
region with input B. All squared coherency curves attained a peak near 
the top end of the input spectrum. Input C consistently produced a 


smaller peak than input B. 
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Table 6.3. Details of the Dynamic Response 


Table 6.3, rows 8,9: and 10, contains estimates of the so-called noise 
residue of the spindle response to the three inputs. These numbers were 
calculated from the respective columns in rows 5,6 and 7 by substracting 
from these the quantities Layo 6) obtained from figure 6.11 converted 
to decibels. The numbers in rows 8,9 and 10 are always greater than, or 
equal to, the numbers in the respective columns of row 4. This indicates 
that the portion of the spindle output that was totally incoherent with 
the input was always greater than, or equal to, the steady-state output. 
This is consistent with the earlier idea that the steady state spectrum 
was involved in the transmission of information concerning the dynamic 
activity of the muscle only to the extent that it acted as a carrier. 
Figures 6.12 and 6.13 show estimates of the gain and phase functions 
corresponding to the squared coherency estimates of figure 6.11. In 
light of the inverse relationship between the width of the confidence 
interval for both gain and phase estimates and the squared coherency 
spectrum (see equations 3.51 and 3.52) the results of figure 6,11 predict 
that the most stable estimates for gain and phase would be obtained near 
the cut off frequency of each of the inputs. Examination of figure Gre 
reveals that increasing the muscle length resulted in an increase in the 


coherent gain of the system. Also, coherent gain in the region above 
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about 1 cps increased at the rate of about 20db per decade increase in 
frequency. Figure 6.13 indicates that increasing the muscle length 
produced a decrease in the relative phase lead of the receptor output over 
its input. At frequencies higher than about 10 cps the estimates of 
the phase function tended to return toward zero, predicting perhaps, that 
the gain function did not continue its upward trend indefinitely. 

At each of the three lengths illustrated by figures 6.12 and 
S13, a transrer function was fitted by eye to the composite picture 
consisting of the three gain function estimates and the three phase 
function estimates corresponding to the three inputs. The transfer 


function was of ‘the form 
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(6.1) 
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f3 was required to account for the +20db per decade slope of the 
magnitude characteristic in the region above 1 cps and fy, was required 
so that the phase function would return to zero at high frequencies. 

f,; amd f> produce the non-zero phase lead present in the system at the 
lowest frequencies. The dotted lines in figures 6.12 and 6.13 are the 
fitted functions. The four break frequencies are sufficient to account 
for the main features of all the estimates. By systematically changing 
f1, fo, £3 and f4, the required reduction in the phase advance produced 
by the spindle with increasing muscle length could be achieved. No 
difficulties were encountered in fitting the same four break frequencies 
to the estimates of gain and phase at a particular length with any of 
the three inputs. It was, however, impossible to maintain the same gain 
(H) and still obtain, a-satisfactory fit. The break frequencies and the | 
gain at .35 cps for the dotted curves in figures 6.12 and 6,13 are 


shown in table 6.4. 
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Table 6.4. Parameters of the Linear Model. 


The coherent gain of the muscle spindle was about 8 times greater (17- 
18db) at the longest length than it was at the shortest length. To 
enable comparison, the curves fitted to the estimates of gain and phase in 
Pisures, 6.12 and, 6.53 input C, are consolidated in figure 6.14. 

Figure 6.15 shows estimates of the squared coherency spectrum 
and the gain and phase functions of a muscle spindle (L = 22mm) under 
two amplitudes of dynamic stimulus. Estimates in the left-hand column 
were obtained using input C at about 45ym rms, while those in the right- 
hand column were obtained using input C at 150 pm rms. Equation 6.1 
was fitted by eye to the left-hand column and the same curve was drawn 
in the right-hand column for comparison with the estimates. Minor 
changes do occur in all three types of estimates. Both the coherency 
and the gain function are shifted upwards slightly by the increase in 
stimulus amplitude. The tendancy of the phase function to return toward 
zero at high frequencies and low stimulus amplitudes is no longer present 
under the increased input. 

Figure 6.16 shows the response of the spindle described in 
table 6.4 at Ltt. Omm to a ramp increase in length. Using the para- 
meters from table 6.4, the response to the ramp predicted by equation 
6.1 can be shown to be proportional to 
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LINEAR MODELS DESCRIBING MUSCLE SPINDLE DYNAMICS. 
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Figure 6.16 


INSTANTANEOUS FREQUENCY OF THE SPINDLE RESPONSES TO A RAMP INCREASE IN 


MUSCLE LENGTH. 
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where t; is the time during which the ramp is applied. If as in figure 


6,16, ty=3 seconds then equation 6.3 becomes approximately 
ae ale (6.4) 


for all t slightly greater than ty. Figure 6.16 definitely confirms the 
presence of a process with a time constant of about 3 seconds. However, 
another process with a much longer time constant is also indicated. 

Data segments for the spectral analysis were apparently too short to 
allow detection of this process. Equation 6.1 can be easily modified 

so that the additional feature of the ramp response made apparent by 
figure 6.16 is taken in account. Assuming that the additional time 


constant is about 50 seconds, the ramp response predicted by 


H@lt526) Gietlos) (14+. 12s) 
(1+50s) (143s) (1+. 002s) 


(6.5) 


will closely follow that in figure 6.16. The ramp response predicted 


by equation 6.5 with t,=3 seconds is approximately 
auger su SW eyes) (6.6) 


for all t slightly greater than t;. At t=t, equation 6.6 shows that 
the amplitude of the 3 second time constant process is about 5.7 times 
greater than that of the 50 second process consistent with figure 6.16. 
Equations 6.5 and 6.6 predict that the spindle will have a small static 
response at the muscle length at which they apply. This is consistent 
with table 6.1. 

The power spectrum of a point process with nonindependent 
intervals can be calculated from equation 5.14 if the characteristic 
functions of all order intervals are known. If the density functions of 
all order intervals are bounded, then all of the characteristic functions 
will tend to zero as the frequency becomes large (96). Therefore, if 


the probability density functions of all order intervals are bounded 
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and the infinite summation in equation 5.14 is finite for some frequency 


fy then for fot. equation 5.14 will tend to 
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Figure 6.2 implies that the density functions of all order intervals 
between muscle spindle discharges were bounded. Indicattons from 
figure 6.10 and from all muscle spindles examined under static stimu- 
lation were that the power spectrum of the sensory discharges was near 
its asymptotic value in the region above about 10 cps. If it is 
assumed that the convergence was to that predicted by equation 6.7, the 
absolute gain of the spindle receptor can be calculated. 

According to equation 6.7 the power spectrum of the sensory 
discharges approaches a quantity which involves the area of a single nerve 
action potential. The shape of the action potential was calculated 
previously (99), and the results of this work were used to calculate A. 
For a 15 pm nerve in Xenopus laevis the area of a single action potential 
at 20°C is about paue y volt-seconds. In figure 6.10 row 3 column 1 
(r3,cl), €=60 msec. Therefore, the power density of the spike train 
should approach about Bron volt?-seconds. This is assumed to correspond 
towabout Wodb ain ficsure 6.10. Hance, at eps tin ficure, 6.10 (3 yen), 
the power density drops to about ome volts?-seconds. With the 
application of input A in figure 6.10 (r3,c2) the power density at 1 cps 
increases to about Cyaan volt?-seconds. Therefore, input A at a 
muscle length of L ti. Onm increased the power density in the spindle 
response at 1 cps by eeerannne volt2-seconds. Input A was measured at 
about 45um rms distributed as shown in figure 6.9A. Therefore, the 
power density of input A at 1 cps was about ie On neeeee ends: Thus, the 
overall power gain of the receptor at 1 cps with input A was about 
90 volts?/m2. Using figure 6.11 (f3,cl), the coherent gain of the 
spindle at the same frequency and with the same input was about 9 volts/m. 
This value was reduced somewhat for inputs B and C as indicated in table 


6.4. Finally, using equation 6.5 and assuming input A, H was calculated 


to be about 1.3 volts/m. 
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CHAPTER 7 


' DISCUSSION 


The results of the previous chapter indicate that methods of 
spectral estimation can be meaningfully applied to the analysis of a 
sensory receptor. By filtering, most problems associated with the 
detection of the point sensory discharges by continuous analog to digital 
conversion were eliminated. The transformation of the discrete sensory 
events into a continuous waveform enabled the analysis of both the stimulus 
and the receptor response with equal facility. Equi-spaced samples were 
obtained from both the stimulus and the response and spectral estimates 
were calculated using modern techniques (80). 

Watanabe (100) has used techniques of spectral estimation to 
obtain estimates of transfer functions between simultaneous nerve impulse 
trains in the crayfish brain. Although it was claimed that the treat- 
ment of nerve discharges was analogous to the treatment of point processes 
by Bartlett (101), this was not the case. Watanabe sampled the spike 
trains by counting the number of discharges in the sampling interval. 

This procedure is equivalent to applying a digital difference filter to 
N(t) the cumulative number of impulses up to time t. Mathematically, 

N(t) is related to the point sensory process by integration. It can be 
shown that the gain and phase functions introduced by the integrate and 


difference operations are respectively 


sin? wea 1 
= ——. sll 
G(f) TE SF, (7.1) 
and 
PCR ee eR pee (7.2) 
a paN 


where A is the sampling interval. Therefore, the power spectrum of a 


point process is apparently linearly related to the power spectrum of its 
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average frequency and exactly equivalent at f=0. All of Watanabe's 
spectral estimates were weighted by equations 7.1 and 7.2. In addition, 
since all the spectra were weighted by the same function, its effects 
would have cancelled in estimates of coherency, gain and phase, Therefore, 
it might seem that Watanabe's transfer functions did actually relate the 
power spectra of the point processes. However, consideration of the 
aliasing properties of N(t) indicates that this is not always so. 

Suppose, as evidenced by discharges from the amphibian muscle spindle, 
that the point process from which N(t) was derived had a flat power 
spectrum in the region around the folding frequency of 50 cps. Then, 
using arguments similar to those in chapter 5, the aliased power in the 
sampled version of N(t) at 40 cps will be about 40%, at 25 cps about 15% 
enaratelOtcps abouty22. oIn general, as stated in chapter 5, it is not 
possible to remove from a sampled signal that portion of the power 

which is aliased from higher frequencies. Therefore, it is impossible to 
relate Watanabe's transfer functions to those which relate the power 
spectra of the point processes except possibly at frequencies low compared 
to the folding frequency. At frequencies where aliasing is negligible 

the two transfer functions are identical. 

If, as in this thesis, the point process is suitably low-pass 
filtered to control aliasing, then the methods of cross-spectral analysis 
described are extremely well suited for obtaining estimates of the 
transfer functions and the squared coherencies between two or more 
simultaneous nerve impulse trains. If all the spike processes involved 
are modified using filters with identical low-pass characteristics, then, 
as pointed out in the previous paragraph and by others (86), the effects 
of the filters cancel from estimates of the transfer functions and 
coherencies. 

Sampling the ayerage frequency of a spike train results ina 
loss of the fine details in the impulse pattern. That is, any information 
that is encoded in the impulse interval or in the relationship between 
adjacent or near adjacent intervals is lost, Also, relating the average 
frequency of the sensory spike train to the stimulus implies that the 


decoding mechanism responsible for the recovery of information from the 
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receptor is concerned with just the ayerage frequency of the spike train. 
Evidence exists which indicates that sensory decoding more closely 
resembles low-pass filtering (13,16-20), a process which is strictly 

not equivalent to averaging. The results of this thesis show directly 
what information can be recovered from the receptor response by various 
filtering operations. 

Comparison of the power spectra of the sensory discharges before 
and after the application of a dynamic stimulus has shown that although 
its spectral distribution was altered, the stimulus was present in the 
discrete discharges of the receptor output at the appropriate frequencies. 
The response to the dynamic stimulus was apparently an addition to the 
response to the static stimulus upon which the dynamic stimulus was 
superimposed. The power added to the receptor response by the dynamic 
stimulus was not entirely coherent with the input (see table 6.3). 
Transmission of dynamic information by the spindle was therefore a 
"noisy' process. In the band of frequencies studied, the portion of the 
dynamic response that was linearly correlated with the dynamic input was 
related to both the amplitude and the velocity of the input. At short 
muscle lengths, the response was almost entirely proportional to the 
velocity of the stimulus. At longer lengths, the spindle became sensitive 
to the amplitude of the stimulus as well. The coherent dynamic 
characteristics of the receptor were only slightly altered by a three- 
fold increase in the strength of the stimulus. Therefore, in that sense, 
the receptor was studied in a linear range of its operation. 

The response of the amphibian spindle to a ramp increase in 
length closely resembled that obtained by Poppele and Bowman (9) from 
the primary ending of the mammalian muscle spindle. In the de-efferented 
state,neither spindle showed much sensitivity to a static stretch and 
both showed a two-time constant decay from a ramp stimulus. The first of 
these was of the order of a few seconds, while the second was considerably 
longer. However, in the mammalian primary ending these time constants 
were not changed by different levels of static length stimulation. 
Changing the muscle length resulted only in gain changes in the mammalian 


spindle with no changes in the relative phase function between input 
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and output. These differences between the amphibian muscle spindle and 
the mammalian muscle spindle primary ending are apparently not the result 
of different analysis techniques. Although the mammalian spindle was 
examined by sampling the average frequency of the sensory discharges, 

the sampling rate was always adjusted to be much greater than the stimulus 
frequency (9). Changes in the time constant of the generator potential 

in the frog spindle from a ramp stimulus has also been reported by 
Ottoson et al.(102). Poppele and Bowman (9) report that the primary 
ending of the mammalian spindle tended to phase lock to sinusoidal stimuli. 
No related phenomenum was observed in the amphibian spindle response to 
random stimulation. The tendency of a neuron to phase lock to a sinu- 
soidal stimulus is apparently related to the variability in the discharge 
before stimulation (103). The more variable the discharge the less the 
tendancy to phase lock. Stein and Matthews (104) report a coefficient of 
variation of about .05 from primary muscle spindle endings in the cat. 
This is almost 10 times less than the variability obtained for the frog 
Spindle discharges. Therefore, although the use of a random stimulus 
probably reduced the tendency of the amphibian spindle to exhibit phase- 
locking properties, it cannot be considered the sole cause. 

It is now evident that the estimates of the gain and phase 
functions of the amphibian muscle spindle could have been improved by. 
shaping the input spectrum. In light of the form of the estimated gain 
function, application of a flat-spectrum stimulus resulted in the highest 
frequency components of the stimulus being the strongest in the 
receptor response. As a result, the coherency between the stimulus and 
response was greatest in that region of the spectrum. Therefore, estimates 
of gain and phase were most stable near the high frequency end of the 
stimulus. If on the other hand, the inputs were shaped to compliment 
the characteristics of the spindle, uniform stability of the gain and 
phase estimates could have been achieved. In this sense, a better 
stimulus would have been one which decreased at the rate of 20db per 
decade increase in frequency. Using the shaped input, an estimate of the 
“best transfer function of the receptor could probably have been 


obtained from the single application of a single input. Also, shaping 
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the stimulus may eliminate or at least reduce the frequency dependent 
nonlinearity present in the receptor. 

Stein and French (103) have shown that the effect of variability 
in their neural analog is to prevent noise-free recovery of encoded 
dynamic information by low-pass filtering no matter how widely the 
dynamic stimulus is separated from the mean discharge rate. This is 
because a random discharge contains power at all frequencies. If this 
random discharge is acting as a carrier of dynamic information, then the 
recovery of that information by low-pass filtering wili at all frequencies 
be contaminated by the incoherent carrier. As the variability in the 
discharges of the neural analog was reduced it was shown that the power 
spectrum of the discharges contained less power at low frequencies 
permitting a more uncontaminated recovery of the information at those 
frequencies. Ultimately, for a completely regular carrier, Bayly (15) 
has shown that distortion-free recovery of the signal is theoretically 
possible by low-pass filtering if both the stimulus frequency to carrier 
frequency ratio and the depth of modulation produced by the stimulus 
are sufficiently restricted. 

The low-noise recovery of the dynamic stimulus from the amphibian 
muscle spindle is facilitated by the fact that the power spectrum of the 
discharges is low at low frequencies. Therefore, information concerning 
the dynamic state of the muscle which exists at these low frequencies 
can be recovered with less contaminative noise than information of the 
same magnitude at higher frequencies. The reduction in low-frequency 
power can be predicted from the distribution of the interspike intervals 
and is not altered significantly by a weak serial correlation between 
adjacent intervals. Coherent information is recoverable from the spindle 
by low-pass filtering at frequencies well past the mean discharge rate 
of the carrier. However a larger percentage of the total signal recovered 
will be incoherent with the input as the frequency is increased. Bayly 
(15) suggests that parallel channels each with random carriers could 
result A an overall reduction of the noise distortion since the carriers 
will likely be mutually incoherent. 


Estimates of the power spectrum suggest that the steady-state 
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discharges from the amphibian muscle spindle can be adequately described 
as a renewal process. That is, in the steady-state, the distribution 
of the first-order intervals constitutes a complete description of the 
point sensory process. Theoretical considerations of the power spectra 
of independent point processes have also shown that the exact form of 
the distribution chosen to describe the renewal process is not critical. 
However, the Gamma density must be preferred over the Gaussian density 
for the amphibian spindle. As shown by Buller et al. (54), extensive 
changes do occur in*the distribution of the interspike intervals from the 
frog spindle at very low muscle lengths. At lengths which result in 
greater than ‘critical depolarization’ the distribution of intervals is 
nearly symmetrical. At lengths which result. in less than 'critical 
depolarization’ the distribution tends to be exponential. The Gamma 
density can be both exponential and highly symmetrical in form. 

The overall dynamic characteristics of the amphibian muscle 
Spindle are usually considered to be the result of the contributions 
from three functional subsystems. The first of these is referred to as 
the mechanical process. Any modification of the stimulus by this sub- 
system is due to the visco-elastic properties of the intrafusal muscie 
fibers upon which the sensory endings reside. The second is the trans- 
duction process. Here, the output of the mechanical process is converted 
to electrical energy as a generator current. Finally there is the 
encoding mechanism. Here, the resulting generator potential at the site 
of impulse initiation is encoded as a series of 'all-or-none' action 
potentials. Vallbo (105) has obtained evidence that in the frog, the 
region of the sensory axon in which the encoding mechanism resides has 
different properties from other portions of the axon. He noted that the 
frequency of discharge elicited by mechanical stimulation of muscle 
spindles was of longer duration and of a wider frequency range than the 
repetitive firing due to electrical stimulation of the sensory axon. 
The evidence suggests that the gross mechanical properties of the intra-~ 
fusal fibers in the frog spindles do not contribute to the dynamic 
characteristics of the receptor except possibly to attenuate the stimulus 
O27, 23,25). tise, the prequency of the action potentials produced by 


the encoding mechanism is apparently linearly and instantaneously related 
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to the magnitude of the generator potential (41,56). Therefore, it would 
seem that the most likely origin for the dynamic characteristics of the 
spindle must be in the transduction process, On the other hand, the 
mechanical process could be involved if the dynamic response originates 
in the microscopic connections between the muscle fiber and the nerve 
terminals or in the mechanical properties of the nerve terminal membrane 
or even in the collective mechanical properties of the terminal nerve 
bulbs and interconnecting cylinders. Certainly, if the dynamic behavior 
is mechanical in origin then, for example, a nonlinear stress-strain 
relationship could account for the reduction in the relative phase 
advance produced by the receptor with increasing levels of static stretch. 

The observation that the gross mechanical properties of the 
amphibian intrafusal muscle fibers do not contribute to the overall 
dynamic response of the spindle does not necessarily imply a similar 
behavior of the mammalian spindle. Indeed, there is some evidence that 
gross mechanical factors are involved in the mammalian spindle (24). 
Also, the finding that different levels of static stretch do not change 
the relative phase advance in the response of the de-efferented mammalian 
primary ending to dynamic stimulation (9) while increased static stretch 
does reduce phase advance in a similar frog preparation, suggests that 
the two receptors function by somewhat different mechanisms. Although 
the secondary ending of the mammalian spindles does show a changing phase 
advance with stretch, this ending is also sensitive to different levels 
of static stretch (9). When the changing time constants are scaled with 
respect to the mean discharge rate, they once again become invariant. A 
similar scaling is not possible with the amphibian spindle since there is 
very little change in the discharge aa tfeay increasing static stretch. 
Invariance of the time constants would not be achieved by scaling with 
respect to muscle length since the time constants do not all change in the 
same direction (see table 6.4). Therefore, the mechanisms underlying 
amphibian spindle behavior must also be considered as being somehwat 
different from the mammalian secondary ending. 

While stretching a de-~efferented muscle it was observed (25) 


that the length of the reticular region of an intrafusal fiber increased 
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nonlinearly. At short oyerall muscle lengths, the strain in the ‘reticular 
zone’ was a very small percentage of the strain in the overall muscle. 
When the muscle was stretched past Ly the length of the ‘reticular 
zone' increased dramatically and then saturated. Further increase in 
muscle length produced no further strain in the ‘reticular zone’. At 
lengths greater than Ly the overall sensitivity of the spindle to dynamic 
stimulation also increased suddenly. It was also noted (25), that when 
the intrafusal fiber was made to twitch under isometric conditions by 
electrical stimulation of the polar regions, the length of the ‘reticular 
zone’ increased during the twitch to its saturation length. Presumably, 
continuous fusimotor activity would therefore maintain the ‘reticular 
zone‘ at the saturated length. This implies that by stretching the muscle 
well past Ly? the effect of motor nerve activity on the sensory ending 
was at least partially simulated. The adequancy of the simulation depends 
on whether the numerous nerve contacts on the ‘compact zones' adjacent 
to the ‘reticular zone' (45) serve simply to anchor the ending so that 
the length of the reticular region can be sensed. In this case, a good 
simulation of intact fusimotor activity would have been obtained. 
However, if the contacts on the ‘compact zones' are sensitive to strain 
in that region, the simulation would have been only partially successful. 
Verveen and Derksen (106) have been able to measure the noise 
fluctuations in the membrane potential of a frog nerve. For a 4um nerve, 
the wide-band rms intensity of the noise was found to be about 200nV. 
This is almost an order of magnitude greater than that expected from 
theoretical considerations concerning thermal noise in the membrane. If 
the wide-band rms intensity maintains this relationship to the expected 
thermal noise then, in the fine .lym nerve terminals of the frog muscle 
spindle, the noise intensity would be about 3mY (the figure for expected 
thermal noise obtained from (54)). This is a significant fraction of the 
15-20mV required for threshold depolarization and could thus account for 
the variability seen in the discharges from the amphibian muscle spindle. 
Presumably, muscle spindles are involved in an animal's ability 
to control and maintain its posture and in its ability to execute precise 


predetermined movements. Certainly, control of posture and movement by 
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Simple activation of the motor nerves with no sensory feedback of the effects 
would be grossly inadequate. Such an ‘open loop' system would not permit 
any compensation in effort for say muscle fatigue or for misjudgement of the 
magnitude of the load. In mammals, it is known that stretching sas fully 
inneryated muscle usually results in a reflex contraction which tends to 
resist the extension. In humans, this ‘stretch-reflext is easily demon- 
strated as the tendon jerk. The reflex is mediated by a sensory volley, 
initiated by the stretch, from the muscle spindles. The ‘stretch-reflex' 
therefore has the properties of a length seryomechanism, that is, of a 
closed-loop system in which negative feedback is used to maintain a constant 
length in the midst of external disturbances. It has been suggested (107) 
that in the mammal, precise movements are initiated by central activation of 
the small motor nerves which innervate the intrafusal muscle fibers (y 
motoneurons). Activity of the y motoneurons excites the sensory endings of 
the muscle spindles which in turn monosynaptically excite the large motor 
nerves which innervate the extrafusal fibers (a motoneurons). This scheme. of 
control is known as the ‘follow-up servo‘ hypothesis (107). The control 
inputs to the system are the y motoneurons. Error between the desired 

and actual limb position for example, is measured by the muscle spindles 

and the error signal is fed to the central nervous system along the 

sensory nerves. The errors are continuously processed by the central 
nervous system and corrective signals are sent peripherally along the a 
motoneurons. The extrafusal muscle fibers act as the actuators which 
correct the error. The ‘follow-up servo’ hypothesis is a very attractive 
idea since any conflict between the voluntary systems and the postural or 
reflex systems which resist the change is simply and elegantly resolved. 
However as has been pointed out (8) there are difficulties with this idea. 
First of all, large modulations in the discharge frequency of the y motor 
nerve produce relatively modest changes in the sensory discharge from the 
spindle, Second, it is not always true that the activity in the y fibers 
leads the activity in the 0 motor fibers. In fact, rapid movements are 
sometimes implemented exclusively by a motoneurons. In light of these 
difficulties, the current idea is that Born o and y motoneurons are 


intimately involved in the control of movement. It is postulated that 
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three distinct controlling signals from higher centers are involved in 
voluntary control (108). The centers separately control the a moto- 
neurons and the y motoneurons. A third signal controls the sensitivity 
of the monosynaptic reflex within the spinal cord. The sensitivity of 
the agonistic muscle system is turned on during voluntary movement 
while the antagonistic system fis desensitized. Deliberate movements are 
controlled by parallel y and a activation and a sensitized reflex loop. 
Activation of the y nerves constitutes a control input and is available 
if compensation against unexpected loading conditions is required. This 
scheme is considered a ‘servo-assisted' method of producing movement (8). 
The 'stretch-reflex' which is apparently of prime importance in 
mammals can only be weakly elicited in amphibia (109). The postural 
activity in the frog appears to depend much more on afferent impulses of 
cutaneous origin. However, a weak stretch reflex does not necessarily 
imply that muscle spindles in the frog do not play a part in the precise 
control of movement. Assuming that the current thinking on the control 
of movement in mammals is correct, then it is not unreasonable to assume 
that the control of movement in amphibia is a simplified version of the 
same thing. Intrafusal muscle fibers in the frog are co-innervated with 
extrafusal fibers. Separate control of the sensitivity of the frog 
muscle spindle is therefore not possible. However, as mentioned earlier, 
it was observed (25) that any fusimotor activity seems to stretch the 
"reticular zone' of the spindle into saturation, Hence, it may be that 
the amphibian muscle spindle operates as a two-state device. If this 
were so, then sensitivity of the sensory ending would not depend on 
the magnitude of the activity in the motor nerves but only on whether or 
not there was-any activity ati all-* Instead of three central controlling 
pathways as in the mammal, it seems that the frog would require only one, 
Activation of the agonistic motor system would fully sensitize the 'seryo- 
assist' loop as well as initiating the movement. Inhibition of the 
Bene tie motor system would densensitize any reflex loop which might 


resist stretch in the antagonistic muscles. 
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APPENDIX 1 


THE DIFFERENTIAL PULLER 
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Figure Al.1l 


CIRCUIT DIAGRAM OF ONE SIDE OF THE DIFFERENTIAL PULLER 


Component List--Ficure Al.l 


Al 50 watt Operational Amplifier 
A2 Motorola MC1439G Operational Amplifier 
HP24DCDT-050 Hewlett-Packard Length Transducer. 


Maximum Displacement +.1 inch 
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Figure Al.2 


CIRCUIT DIAGRAM OF THE SIGNAL SPLITTER. 


Component List--Figure Al.2 


A3} 
A4 


Motorola MC1439G Operational Amplifiers. 
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Figure Al.3 shows the gain and phase performance of the 
differential puller. The gain and phase relations are shown between the 
control input and the output of the HP24DCDT-050. The dotted curves show 
the predicted relations between the control input and the actual dis- 
placement produced by the puller. These curves were obtained from the 
solid curves by applying the correction formulas for the HP24DCDT-050 
supplied by Hewlett~Parkard. 
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Figure Al.3 


FREQUENCY RESPONSE OF THE DIFFERENTIAL PULLER. 
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APPENDIX 2 . 


THE PULSE-HETGHT DISCRIMINATOR 


Component List--Figure A2.1 


A5 Motorola MC1439G Operational Amplifier 
A6 Minimum pulse amplitude } 


A7 Maximum pulse amplitude | Fairchild p»A710C Comparators 


G1-G10 Motorola MC724P Nor Gates 
FF1-FF2 Motorola MC790P j-k Flip-Flops 
BL Motorola MC/799P Buffer 


The output of Bl falls from the 1 state to the 0 state 
coincident with the fall of the input pulse if and only if comparator A6 
was set by the rising phase of the input pulse and A7 was not. 
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Figure A2.1 


CIRCUIT DIAGRAM OF THE PULSE-HEIGHT DISCRIMINATOR. 
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APPENDIX 3 


THE SPIKE FILTER 


Pieure AS. k 


CIRCUIT DIAGRAM OF ONE SECTION OF THE SPIKE FILTER. 


the transfer function of the circuit in fisure A3;) is 


-R ; 
ee ee (A3.1) 
al ee ee 5 
Se yn ene 
1 + [C)RoR3 e Ro R/S 2R3C1Cos 


If the components are chosen as shown in figure A3.1 equation A3.1 becomes 


| ¥2s , s@ (AAS 1) 


The characteristics given by equations 5.27 and 5.28 were obtained 
by cascading two stages of the circuit shown in figure A3.1l. 
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APPENDIX 4 


THE CROSS-SPECTRAL ANALYSIS--FORTRAN IV SOURCE PROGRAM 


DIMENSTON X12048),4(20481,Z{ 802,4),F0 802),F AX (36), BATA( 2048), 
Hee (22) Ale BOOB BO2),IS11024,2) yAVEL29 sTITLESS 5 2) NEES IS 
ee 6) VARIO EI eG 5) oy ae iY Ao 
DATA XOI™,YDIM, OPTS, FPEQM/4.0,4.050.07;23.0/ 
DME A MOTTE G/4C34ERENCE RAS sue ceaEs MAGNITUDE — 
eee ee ee te 
NATA havo, 155147 
CAlS PL aT SUBATAL 1398102) 
NEEC Sey, is oe 7 a = a 
Be 28a eas 07 
eee alee es 
lial = ON He ey Oa, 
25 CONTENUE 
Gm, 20 = 1,4 
O21 d=1, 802 _ a bee Tue oF ; ‘ 
TEAO PR lk eee 
21 CONTINUE, 
a0 CONTINUE 
PCR PAP G(G) 
Pees 1 Naess 
HEC O UNC vels UCL OS rile el es aces ee Lae eee 5 Bye 


5 Ie pe: 


Weis Gtr. 
READ (1,190} 
BWET LEE. 35) 
Ree uel dee 
Diy 22 eee 
WB CONTINUE 
et ee Oee ee eri 


Shee 


| esa a 


55 
AVE J=0. 
VAR(JI=0. 
NO 12 1=1,5 
eo = et 5) 
ale ly JY 
XCLtS1L 21274 
See eo iehs 

PO 10 J=1,1 
BE) Se om 


Aae 


AVGTJV=AVS TS 


Y(I)=0. 
10): CONT INGE 
DR T=192 


20 


wGLPSO. 
YUL P20 


Or Ai Ea ese 


=FLOAT(S 


Pe Fo dl, 32) 
GA G12) 5 eee lad a 


= a. 
Mie Ula )ia Ja ts 2 9 2 


2) 


71 


? 
deacl 


ne 


(JS EE ASS Wee 


Los) i 


tS 


Es we er PON da ee ee ore amed ee 
Oise oS 00 eo), l= 140 812) 


Bled 


tiple 


ies 
P50 


O24 
in hae Go ie] 
Nate Bey {A rete mn sere Weis 


5 


2048 


ye 


aye 


UGA GI Pinos atl ef 
UE 30 tee Baa ee OL 


AC OR tM Gales leat COURS S| 


Petia eet baa ty 4 Se tle ta tS a es ae 7 a — 
B= OWS elle COS (Pree Opt Peb 37 10234-7) 
rl De dec a fd OR 
OC CaN 
VARI J)VSVAR (7/1023. 
eS PSA ON 14 e x , Vg 2) 
Ca oi t= 1 80? . | 
MEL yeaxt 1) 
ee eet in) 
Go Cen TIVE 
550 CONTINUE 
Die cc be Ea ey 00 Baek Oh Bn ic! eae) m 2. te 
ee ty a) (NOL) it, 2) 
Dehes- T=, 80 
Cl EvAL Rees Vy OGhe set heh oan 
Pe et PN YN EE aA Tp oRY tT) 
CE) a a Pe Cd eh EK 8 ET PRM OTD) 
PAA yee CLS yt ( 1 pe deme: ME eae = een Ae 7 
&3 CONTINUE 
Sel (MENA LEME 
ee i at ee 
PQ 2 f=] 4% 
pee ep a) Se UDA ra ol) PO. ss t= 3) Ae ve as. 
3 CONTINUE 
2° CONTINUE 
ied 
“8 eee ee 
[FL FAXTKK)..LE.FREQM)GO TO 5 
Nihics ee eas et nd oe A SENSED ees . br 
Po emote Lie See) 
FSPH=LALOGLOL FO) —FSI9ISXDIM 
1=2 
IT=0 
eee a Ue ot 
Bee ya es Op), eel) aes i ee ee be od en 1 
ES Se Sa 
34-FSFeTALNGIOIR(J2)-FS1I/FS2 
Dosti Sih iat aio ert ia te oa) 
EST Sece 
ed oe 
Tetris. eee Ge PG 35. a an = 
GC TR 34 
eet [pod 
Leet Oe om lt. oe 
NEOr Riga fji—h)/4 el 
PCISI=ACESTHRSIT S23 
I=J-] 2 = _ 


Woe, Saari A ee or 
ap Te ! aur tea ap 
pe ‘ ie) pte a — nines 
a: os . ; ; 
| , 7 ee | }3u= 
rn ea a yt eb Trt saan ay 
bE Por A de pp prerNte e Bens | 
cs . AMG Ya - | 
mi eae a tt FS 
Rea am uth heave tui Ls, 
AM iy A Pa 
te , teeta ner bik 
| ! é as rane i an) of in 
i ey Tie > Ale 
— ; > { iyvat & i 
‘ a Fe > do ny ays Fe 10 
: ‘: | span aes 
a ee oe ae - ite a Ki 
re fd PO ANOS Ea att Ry’ 
; S bok ae ae as 
i* re th SAmes ot i 
RR es SER SL) | 
ik ten shade, ieee hy 
| OTe A ee nght felons 1 
3 HER c o METS ree ne 
| 1 Caer 
“ln 
; Os ial) S a 
een te £7. 
si Lee rs | ; eng’ i. ike 
h - 
- 1 to 
A 
{ he 
‘ ey LF 20-4 
Tae 
~ 4 en _ 
ee ae semanas ia 


ely “A -. ihe ma bret 


ie ee ¢ ba Aree 2 petnies 
iY rey a con iy ap 43% 


5 


eC Ee vaee Uae 
Z1=Z14+7(L,™) 

CONTINUE 

(TI MIEZIUS FLOAT ISH1L 1) L 


CONTINUE 
f=J4+1 


BOUPe UL aera 


DB 65 


Sea eth 
Sh = Sie ee alae Dy ek 4? 


aS. © fo eee XXe7 


id 
‘iu 


12 VRKZ) 
AV) 


Ee a9 ] ee 
Zid y 2 oRTAN2 EFL 1,3).2 pe AC? wee ee UCol y tae) 
Po ys ae Pe uty a Cl ye IP | 
INA eo? Pa se Vee, (PT ee ewes ae 2 ; ee: 
71143) = SKY /SXN% (Ch IAC T EE) 

5S. CONTINUE 
Ue ie a) eS 8 oat 
GALT WY 2UGs (C4057 S05 =3) 
Vile oe 
etek ae eee : : at, jo 5 bes We sa m! rine. 
el ite2iat) 
T=KK-G 

ne een a 5 

Lek: Whee AGT 
TPA) RE. E/ RAR AGE NT e 47 
KR=KA—-4 st oe , SSR ae 
eer ] 4, es Ge 3 arf 19 ° 

Ciena is 
AGG Gare se v JYIIM 
Z(114#268,1)=0. 
LAI IAZS Sl SLt 1 14 2.) 
AeehiG wei s ela uct ee aed Pee cou (eet pear? 
LOL IESSe Tle epee yt) 
LEA TEL eI SH 20 
ete Me et ig eC ae, Se ke Ui kee Se Sa : 

GS 2a 2 0 i5 2) 
FSF=FST 

tr en ae ee 4 ee ace ae ze eevee LS ts eee 
BS PaAMTNIAES 1474 be 309) 
PSEHAMAXLUEUPSEsZ (ly 31) 

Peete le ee Oe Oe noe Wer A: 
RS5h=0,59%FS: 
PSee ll, SE Sk 
‘AGRE ral 22) = ip are pars eas —— - 
MER E23) =26 OF ALOGLO(ESFE/ESI)/YDIM- 
nC (652 ein + 

Qui NCES COLO. 0,2 tL eel isle. 2 ie Or OO oe tek I en 

LUT i+ M)520.9 
GAL >LOT (OsRyO 20 ss) 


=159— 


7 
fi 
~ ‘ ‘ ' } 
Pe, 
¥ 
t 
ye 
{ 
— 
ee 
af 
jan 4 
sata > A 
i i 
Me 
Aled ’ ry pe 
aL ! 


aii 


em 


> 


, OE "* ane 


i ih 
eG ‘~, 
x sre 
a i 
a 
Se 
eS kl oe 
ae | ‘4 i 
A ra 
ao 
i 
. be 
} ‘ 
ys? 
‘ ef 
' ee ps 


ia 
Boi. ‘oe . 


Rare vette FOIE 


ea ctian iy 


it Linthe hte 
OREN a 


= 


Kx? 


Bi eG, mt a 
oS Roe oceania CRIN IER 


ONDE hee On? ee 


WE AIG BG KK-ASBND. I.E 0. 140KK)/241)6G0 Th 53 

Meee Bee. PERE TRASK ADCO FO 53 

(i Wye NING Sy 2 0D 0.15 1350.6, -1) - Shit. 

apie dene. (na pe 
503 at Soar eee Pep tcts iteiO BO, — 1) 

ELL .6T.9)K=1 

tke “wyunee (eee 0 ioe. PAO, PAX UL oon 
45 CAM PEt 55,060.03) 

& CONTINUE 

LE CM) ALAS. LS 
Cie om fol T=] ne 

ESPH=EDE*DEDIT) 

Br Mares a rag OD bee ere 

[E( DAC GE sd od OA=0.99999 

Deas Set woe (1 tA / (in -DAD) 

Dlg Cog upewi reg uA . oe Sen | (SF) 3 

eed cutee lane iste tote ry Lo) =O 6 

laietar relic get Pte OGY 4 ene 

714534, LISTANTIDA+1.65/SORTLESTI). 

S(O ote 7a gt I ee ieee nae Al Oe lie eo eae Be Ae 

fl et oe, Dee Ce OSA 4y | pace? 
Pe hugee RIMES ony eee OC a OR 7 aie Bele E i a 

CALI Ae Sn aes ar 7(1;, Py ott yby1s 132) 

Pal Vat a ed 16S) i> Soe %,0) 

VEAL. VINE Ch, 2193521) a Foy by 05.09 

GO. FB SG 
CO ta ee do ett 

ee. 4 Zs a re 


Pe ee 7, 1) =O. OO0ONT 
xT) eSDATED A HUFET=20 185.154 Me Reeve Gir ai Wake eer Eb) 


(T)=xX(T) 


Te ee }=0.99999 
DA =ARSIN(X (I) )#2690./PI 

ee igo el ea ee A . i << eee 
Tf 14+ BG OV eT nw. DG 


ae 


74 
a 


CONTINUE 


Teel Aaleg ie) oA. 


Z. 
Vis SO Ae Riles Ss BS ae G E.1LG39 68 Tp oie, 
CALL SYMBOL ( 


eet 2 bala Pelee Ase Os oss 


Kaa 


DE lik ocaoes # - The outs 
{ 157 gE Sm ae re tibia he Ves ee es Bae Vio wco TA acy. 


a 
Phe ee renee. Q\eldl ial 2) 3/201 it2 7233 6) 


TE(Z 
ALL 
4 pececneernewecerccsnerenneeryrensres aa eciie as Od te Oe Bernama a eer ewswsecs pocorn ease mares meen peer oe 


i powacearerersnewncesnenapdc dann american ansararaseere 


Cor (ea 
K=3 


{ 
G 
Z 
\ 


CONTINUE 


K=3 

as ys Pag 1 

Pepa AG ae, Ole ee Ne (lee) Ft 10. GbR Tm 75 
~114-° 


ae DP key “i vial 


f Li eae oe 
*) r 
ot ee, na | at 
‘ i Fay Ww ae aay els Ath be 
: , co. , Pe ’ aoe 
. ie’ v4 he * Sank he 
p i hi ee oy ere’ ep 
: rj { sep £ 
Myr < ‘ ; ue » ; 
‘ bee ‘ i oie Vis 


aes ier 


1 Mak tal MOL 
hgeth EGE Tog ae a 

H -_ x ‘ a a. we ee * Tete, Cie iy 
Siem: ee ff Vat LY 


y et ot ee Set 
Y ’ Rw} te 


: j rr ee ae 
Ms Fae | 


(Se Waa ie Tye Se re 


ae 148; ish 


4 " 
5 4 “ae Ba i 
: mh ie 2 <1 
M ae i 
} } Oar = 
{ ' I \ » i . “n f ® 
a : me Y 
7 i oy pug 
Cp 


se, aia: hus Ae ke 


aoe e ‘ La hg Al AP en Bp tyryy 


es sa ay it ae iy 
25h Birdieoa ss sii etd ie +(e 
1 Pa MPT ha a tia 


pate te 


of oH r ne 
oo - _ - ar y —s : en she 


i fe Ge tt ihe ae aie wieat a5 racy 


ih ced cid vet ipsiities' i 
ana: a 5 anal a ne 


Hans 
ay 


7 


ane 


aoe 
ee) + UB Ue 
75 K=23 7 ale Po a ae Sar re = - Ss nEEEESAIEDEEEERTEREEEEEEEEEeeee a . — £ 
76 CONTINUE 
COI Tite SG 
4h pA ie ae Tat 
Bea GT) Lee Ue 230% 1.-X(L)) 
= puss, Naz ethyl) Hl 2, ee eee ene 
27 GONT TINUE 


DO 32 cae Merit tae | 


ets ee re ah ela): 


Die ner. 2) =76.9 241 801 O17 (1 +267,3)/FS1) 


7 (£4534, 3) =@0.0*AL0G10(2(L4+534%3)/FS1) 
CONTINUE 
GAGS eS ee eel he 2 1 32) 


tae 

Bae Pd 
BE 267. ASL Eo 
CALL PLOTIUETII,(2( 


Pw A Oe UE i yar ee AOS oe CE ee Gay hes 
Ze ee be era G) | eee oh) 


el LE eed Wie ed oc ae 8 es OD a 
T 


2 
Pie a Ge nad tT 


KR 
CONT I MUSE 


ante 
LRG e ee om Ge Bree at! 


gee ere ee LOT P3324, 2 nGhact ide, SS DIM ae aes 
E ipa) are) ea) ot 


a5 


oa9) 


WA 


ey 


390 


ect an t| =) en ne 

TE(M.EQ.39G60 TI 42 

ee =YL i Ee? ieee peeeaee: = 
TFCYL+VYNIM.GT.2°.)G0 TO 55 

CALL eee ee 73) 

aor = VL — vA yy TM ct 

PAG PIR SSR e Yes Sate Gb wean ae 

Den Coen: fel eae bE ol 

CALL PLOTIO.0,9.0,999) 
ECRMAT(15(A4,A2)) 


PECRMATI LHI // 5M ASL AS ACD) 


FORMAT (1616) 


SORMAT L/S ISI I/20 to) 
SQCM ATI S/S / 71 / Ae 
BERMATULIS1/5% TAVEPA 
mls creer 

oe 

END 


Peer Seas S Pe oe ee ee i 
= TMDUT LEVEL? 9220039 5%) (AVERAGE SUTPUT L 


10.4 
: 
(2 = 


=115- 


a ri had Pig We a iene 


ees. 
| ne 
. " eg | Te oh, 
he wertes act 


es m a 7 a 


we: tow) is 
Oe aya Sia peas fee _ 
~ eeu + 
; es sh) 
\ wee (< I 
rt i Se mt Ni ‘.. 
Petia, 
i” ar ¢* , la ah a 
27 Tunas 
s oF : Pl 
; I i 
é : 
tT: 4 —_ rin 
? ~ 
é j ° f > f "Ned : pa 
OT heh E417 PP Tees 
é », 4 
; 2 eninge tod Ors. 
9 a : 
q i a 
: var 
cate 
<< a a ee 
“TV tet ath 
del ea ee EAS 7a ab | 
/ 7 a ‘ tT) : pA Lie A 
i 7 ay - 
a Ae 
e U & 
ipod ae Ys 
, U (=a 8). 


ae CTE, ¢ 
As / p 4 Su i ee 
, 7 ae sf 


a F oe eis ub Made a ie 


a pases ok aa eos 
oie! Ue Pa es . 
28) | PER 


, natin lek ‘3 


ah aa 


i 
1 : 
; Baie " 
ie mit 
‘ 


“s Ph \ if 


: i 
' 
toad Th yea yy 


i 82990! 


